Introduction¶
This code generates the supplemental figures for the paper "Exploring controls on the timing of the phytoplankton bloom in western Baffin Bay, Canadian Arctic".
Directory structure¶
binder, where is the .yml file for the conda virtual environment.
Environment¶
This notebook was run in a conda virtual environment. The following commands will automatically create the environment and launch the Jupyter notebook:
cd <my_path>/gud_groups
conda env create -f binder/environment_timing.yml
conda activate timing
cd utils/python/MITgcmutils
python setup.py bdist_egg
cd build
pip install MITgcmutils
cd ../../../..
jupyter notebook
Once the environment timing is created, the Jupyter notebook can be launched simply with:
cd <my_path>/gud_groups
conda activate timing
jupyter notebook
First, load libraries¶
Uncomment and run one time the cell below if you have the error:
ModuleNotFoundError: No module named 'MITgcmutils'
# import os
# import subprocess
## make sure the current working directory is gud_groups
# wd=%pwd
# os.chdir('utils/python/MITgcmutils')
# subprocess.call(['python', 'setup.py', 'bdist_egg'])
# os.chdir('build')
# subprocess.call(['pip', 'install', 'MITgcmutils'])
# os.chdir(wd)
%matplotlib inline
import cartopy
import cartopy.crs as ccrs
import collections
import datetime
import gsw
import gudinfo
import IPython.display
import math
import matplotlib as mpl
import matplotlib.colors
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import netcdf_tools
import numpy as np
import os.path
import pandas as pd
import plotly.express as px
import plotly.graph_objects as go
import read_mitgcm
import scipy.interpolate
import scipy.io
import statsmodels.nonparametric.smoothers_lowess
import TaylorDiagram
import vstats_pos
Parameters of the figures¶
# dpi='figure' # default dpi
dpi=600
# ext='png'
ext='jpeg'
Select year¶
years=range(1,11)
nbyears=len(years)
iyeartempo=nbyears
first_year=np.arange(0,365)
last_year=np.arange(365*(nbyears-1)+0,365*(nbyears-1)+365)
array1d_iT1y_iT=last_year
first_year366=np.arange(0,366) # first year for heatmaps
Model¶
Model: traits¶
Set indir.
indir_exp0='data/DataS6_output_mitgcm/exp0'
Create RunInfo object.
info=gudinfo.RunInfo(indir_exp0)
Maximum growth rate ($P^C_{max}$, $\mathrm{ d^{-1} }$).
stod=86400 # seconds to day
pcmax=info.f['pcmax'].to_numpy() # s^-1
pcmax=pcmax*stod # d^-1
Cell nutrient uptake half saturation constant on nitrate ($K_{NO_3}$, $\mathrm{ mmol\ N\ m^{-3} }$):
$K_{NO_3}=a_{K_{NO_3}}V^{b_{K_{NO_3}}}$
nplk=info.nplk
K=np.empty(nplk)
K[:]=np.nan
a_kn_no3=info.params['a_kn_no3'] # in mmol N m^3
vol=info.f['vol'] # in um^3
b_kn_no3=info.params['b_kn_no3'] # unitless
for iplk in range(0,nplk):
grp=info.f['grp'][iplk]
if (0<=grp and grp<8):
K[iplk]=a_kn_no3[grp]*vol[iplk]**b_kn_no3[grp]
Cell minimum stoichiometric quotas of nitrogen relative to carbon ($Q^{min}_N$, $\mathrm{ mmol\ N\ (mmol\ C)^{-1} }$)
$Q^{min}_{N}=a_{Q^{min}_{N}}V^{b_{Q^{min}_{N}}}$
nplk=info.nplk
Q_min=np.empty(nplk)
Q_min[:]=np.nan
a_qmin_n=info.params['a_qmin_n'] # in mmol N (mmol C)^-1
vol=info.f['vol'] # in um^3
b_qmin_n=info.params['b_qmin_n'] # unitless
for iplk in range(0,nplk):
grp=info.f['grp'][iplk]
if (0<=grp and grp<8):
Q_min[iplk]=a_qmin_n[grp]*vol[iplk]**b_qmin_n[grp]
Cell nutrient uptake rate of nitrate relative to carbon ($V^{max}_{NO_3}$, $\mathrm{ mmol\ N\ (mmol\ C)^{-1}\ d^{-1} }$)
$V^{max}_{NO_3}=a_{V^{max}_{NO_3}}V^{b_{V^{max}_{NO_3}}}$
nplk=info.nplk
V_max_d=np.empty(nplk)
V_max_d[:]=np.nan
a_vmax_no3_s=info.params['a_vmax_no3'] # in mmol N (mmol C)^-1 s^-1
a_vmax_no3_d=a_vmax_no3_s*86400 # in mmol N (mmol C)^-1 d^-1
vol=info.f['vol'] # in um^3
b_vmax_no3=info.params['b_vmax_no3'] # unitless
for iplk in range(0,nplk):
grp=info.f['grp'][iplk]
if (0<=grp and grp<8):
V_max_d[iplk]=a_vmax_no3_d[grp]*vol[iplk]**b_vmax_no3[grp]
Plot: traits¶
type_id=np.arange(1,52)
pos_Prochlorococcus=np.arange( 0, 1)
pos_Synechococcus =np.arange( 1, 2)
pos_smalleuk =np.arange( 2, 4)
pos_cocco =np.arange( 4, 9)
pos_diazo =np.arange( 9,13)
pos_Trichodesmium =np.arange(13,14)
pos_diatoms =np.arange(14,23)
pos_mixodino =np.arange(23,33)
pos_zoo =np.arange(33,49)
def make_plots(axs):
box=dict(pad=5,alpha=0)
# --- CELL NUTRIENT UPTAKE HALF SATURATION CONSTANT ON NITRATE
ax=axs[0]
h3=ax.plot(type_id[pos_smalleuk]-2, K[pos_smalleuk],'gh',
markersize=12)
h4=ax.plot(type_id[pos_cocco]-2, K[pos_cocco], 'b+')
h7=ax.plot(type_id[pos_diatoms]-2-5, K[pos_diatoms], 'ro')
h8=ax.plot(type_id[pos_mixodino]-2-5,K[pos_mixodino],'md')
h9=ax.plot(type_id[pos_zoo]-2-5, K[pos_zoo], 'k*')
ax.set_xticklabels([])
ax.xaxis.set_major_locator(mticker.MultipleLocator(10))
ax.xaxis.set_minor_locator(mticker.MultipleLocator( 1))
ax.set_xlim([0,43])
ax.set_ylabel('$K_{NO_3}\ (\mathrm{ mmol\ N\ m^{-3} })$')
ax.set_yscale('log')
ax.set_ylim([0.1,20])
ax.grid()
plt.text(0.05,0.92,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
# --- CELL MINIMUM STOICHIOMETRIC QUOTAS OF NITROGEN RELATIVE TO CARBON
ax=axs[1]
h3=ax.plot(type_id[pos_smalleuk]-2, Q_min[pos_smalleuk],'gh',
markersize=12)
h4=ax.plot(type_id[pos_cocco]-2, Q_min[pos_cocco], 'b+')
h7=ax.plot(type_id[pos_diatoms]-2-5, Q_min[pos_diatoms], 'ro')
h8=ax.plot(type_id[pos_mixodino]-2-5,Q_min[pos_mixodino],'md')
h9=ax.plot(type_id[pos_zoo]-2-5, Q_min[pos_zoo], 'k*')
ax.set_xticklabels([])
ax.xaxis.set_major_locator(mticker.MultipleLocator(10))
ax.xaxis.set_minor_locator(mticker.MultipleLocator( 1))
ax.set_xlim([0,43])
ax.set_ylabel('$Q^{min}_N\ (\mathrm{ mmol\ N\ (mmol\ C)^{-1} })$')
ax.set_yscale('log')
ax.set_ylim([0.004,0.1])
ax.grid()
plt.text(0.05,0.92,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
# --- CELL NUTRIENT UPTAKE RATE OF NITRATE RELATIVE TO CARBON
ax=axs[2]
h3=ax.plot(type_id[pos_smalleuk]-2, V_max_d[pos_smalleuk],'gh',
markersize=12,label='picophytoplankton')
h4=ax.plot(type_id[pos_cocco]-2, V_max_d[pos_cocco], 'b+',
label='other nanophytoplankton')
h7=ax.plot(type_id[pos_diatoms]-2-5, V_max_d[pos_diatoms], 'ro',
label='diatoms')
h8=ax.plot(type_id[pos_mixodino]-2-5,V_max_d[pos_mixodino],'md',
label='mixotrophic dinoflagellates')
h9=ax.plot(type_id[pos_zoo]-2-5, V_max_d[pos_zoo], 'k*',
label='microzooplankton')
ax.set_xlabel('Planktonic type ID',fontsize=23)
ax.xaxis.set_major_locator(mticker.MultipleLocator(10))
ax.xaxis.set_minor_locator(mticker.MultipleLocator( 1))
ax.set_xlim([0,43])
ax.set_ylabel('$V^{max}_{NO_3}\ (\mathrm{ mmol\ N\ (mmol\ C)^{-1}\ d^{-1} })$')
ax.set_yscale('log')
ax.set_ylim([0.005,1])
ax.grid()
ax.legend(loc='lower right')
plt.text(0.05,0.92,' c ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
# --- POSITION
axs[0].set_position( [0.10,0.69 ,0.85,0.27])
axs[1].set_position( [0.10,0.39 ,0.85,0.27])
axs[2].set_position( [0.10,0.09 ,0.85,0.27])
with plt.style.context('mplstyles/traits.mplstyle'):
# Plot
fig,axs=plt.subplots(3,1,figsize=(16, 12.32))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_traits.'+ext
plt.savefig(outfile,dpi=dpi)
Supplemental Figure $<$sup_traits$>$. Parameters controlling the half saturation for growth on nitrate ($k_{NO_3}$), $\mathrm{ mmol\ N\ m^{-3} }$, Fig. $<$traits$>$c) for each numerical type (each numerical species). (a) Cell nutrient uptake half saturation constant on nitrate ($K_{NO_3}$, $\mathrm{ mmol\ N\ m^{-3} }$). (b) Cell minimum stoichiometric quotas of nitrogen relative to carbon ($Q_N^{min}$, $\mathrm{ mmol\ N\ (mmol\ C)^{-1} }$). (c) Cell nutrient uptake rate of nitrate relative to carbon ($V_{NO_3}^{max}$, $\mathrm{ mmol\ N\ (mmol\ C)^{-1}\ d^{-1} }$).
Grid¶
drF is the r cell face separation, meaning the thickness of each depth layer (in m).
It corresponds to delR on https://mitgcm.readthedocs.io/en/latest/getting_started/getting_started.html#grid
gridfile='data/DataS6_output_mitgcm/exp0/grid.t001.nc';
drF=netcdf_tools.read_netcdf(gridfile,'drF')
RC is the r coordiante of cell center (in m). We switch sign to have positive depths.
RC=-netcdf_tools.read_netcdf(gridfile,'RC')
RF is the r coordinate of cell interface (in m). We switch sign to have positive depths.
RF=-netcdf_tools.read_netcdf(gridfile,'RF')
Observations¶
Observations: sea ice concentration¶
The values of observed sea ice concentrations were boolean values of 0 and 1 because sea ice melted suddenly in only one day. Here, these values were retrieved from the file of the forcing fields of sea ice concentration. We made sure these values corresponded to the expected values.
infile_sic=\
'data/DataS3_observations_Qikiqtarjuaq/siarea.qik.2016.365.32bits.bin'
array1d_iT_sic=np.fromfile(infile_sic,dtype='>f')
Sea ice break-up was on 18 July 2016, doy (0-based)=199.
obs_breakup=next(i for i,v in list(enumerate(array1d_iT_sic)) if v==0)
assert obs_breakup==199,\
'the day of year (0-based) of break-up from the file %s (obs_breakup=%r)'\
' is different from the day of year of break-up'\
' from Oziel et al., 2019 in Elem. Sci. Anth. (July 18=199).'\
% (infile_sic,obs_breakup)
Observations: snow and ice thickness¶
The file Ice thickness from https://www.seanoe.org/data/00487/59892/ contains $in\ situ$ snow and ice thickness at the Green Edge sea ice camps 2015 and 2016 ($67.48^\circ N$, $-63.79^\circ E$). Variable 'sample_thickness_cm_average' is the snow or ice thickness ($cm$).
Observations: snow and ice thickness: load dataset¶
def load_thickness():
fname='data/DataS3_observations_Qikiqtarjuaq/66407.csv'
df = (
pd.read_csv(
fname,
dtype={
'mission': 'string',
'date': str,
'date_time': str,
'latitude': np.float32,
'longitude': np.float32,
'sample_type': 'category',
'sample_thickness_cm_average': np.float32,
'sample_thickness_cm_sd': np.float32,
'pi': 'category'
},
parse_dates=['date','date_time']
)
)
df['doy']=df['date'].apply(lambda x:x.timetuple().tm_yday)
df['sample_thickness_m_average']=df['sample_thickness_cm_average'] \
.apply(lambda x:x/100)
return df
thickness_df=load_thickness()
Observations: snow and ice thickness: select¶
We select the observations of snow thickness at the ice camp 2016.
snow_thickness_df=thickness_df[(thickness_df.mission=='ice_camp_2016') \
& (thickness_df.sample_type=='snow')]
# snow_thickness_df[['date','sample_thickness_m_average']]
The start of the period without observed snow cover was defined as the start of the melt pond period. It was on 15 June 2016 (Oziel et al., 2019 in Elem. Sci. Anth.).
obs_nosnow=166 # 0-based
We select the observations of ice thickness at the ice camp 2016.
ice_thickness_df=thickness_df[(thickness_df.mission=='ice_camp_2016') \
& (thickness_df.sample_type=='ice')]
Observations: nutrients¶
The file Nutrients from https://www.seanoe.org/data/00487/59892/ contains $in\ situ$ nutrients at the Qikiqtarjuaq sea ice camps 2015 and 2016 ($67.48^\circ N$, $-63.79^\circ E$). Variables 'no3_um_l', 'po4_um_l', and 'sioh4_um_l' are the nitrate, phosphate, and silicic acid concentrations, respectively ($\mathrm{ \mu }$M).
Observations: nutrients: load dataset¶
def load_nutrients():
fname='data/DataS3_observations_Qikiqtarjuaq/66412.csv'
df = (
pd.read_csv(
fname,
dtype={
'mission': 'string',
'date': str,
'date_time': str,
'latitude': np.float32,
'longitude': np.float32,
'sample_type': 'category',
'sample_source': 'category',
'snow_thickness': 'category',
'code_operation': 'category',
'station': 'category',
'station_type': 'category',
'cast': np.float32,
'bottle': np.float32,
'bottom_depth_m': np.float32,
'target_depth_m': 'string',
'depth_m': np.float32,
'filter_type': 'category',
'quality_flag': 'category',
'dilution_factor': np.float32,
'no3_um_l': np.float32,
'no2_um_l': np.float32,
'no2_and_no3': np.float32,
'po4_um_l': np.float32,
'sioh4_um_l': np.float32,
'n_ice_samples_ice_camp_2016': np.float32,
'method': 'category',
'pi': 'category'
},
parse_dates=['date','date_time']
)
.rename(
columns=dict(
depth_m='depth',
no3_um_l='no3',
po4_um_l='po4',
sioh4_um_l='sioh4'
)
)
)
df['doy']=df['date'].apply(lambda x:x.timetuple().tm_yday)
return df
nutrients_df=load_nutrients()
Observations: nutrients: select¶
We select the observations of nitrate at the ice camp 2016 in the water column made with a GFF filter.
nutrients_df=nutrients_df[(nutrients_df.mission=='ice_camp_2016') \
& (nutrients_df.sample_type=='water') \
& (nutrients_df.filter_type=='gff')]
Find duplicates.
g=nutrients_df.groupby('doy')['depth'].value_counts()
g.where(g>1).dropna()
Series([], Name: count, dtype: float64)
No duplicates.
Observations Qikiqtarjuaq 2016: chlorophyll a¶
The file Pigments, nutrients, Chlorophyll a and Phaeopigments (concentration) from https://www.seanoe.org/data/00487/59892/ contains $in\ situ$ chlorophyll a at the Qikiqtarjuaq sea ice camps 2015 and 2016 ($67.48^\circ N$, $-63.79^\circ E$). Variables 'conc_mg_m3' is the chlorophyll a concentration measured by HPLC (mg Chl $\mathrm{ m^{-3} }$).
Observations Qikiqtarjuaq 2016: chlorophyll a: load dataset¶
def load_pigment():
fname='data/DataS3_observations_Qikiqtarjuaq/66417.csv'
df = (
pd.read_csv(
fname,
dtype={
'mission': 'string',
'date': str,
'date_time': str,
'latitude': np.float32,
'longitude': np.float32,
'sample_id': 'string',
'sample_code': 'string',
'sample_type': 'category',
'sample_source': 'category',
'snow_thickness': 'category',
'leg': 'category',
'code_operation': 'category',
'station': 'category',
'station_type': 'category',
'cast': np.float32,
'bottle': np.float32,
'target_depth_m': 'string',
'depth_m': np.float32,
'pigment': 'category',
'conc_mg_m2': np.float32,
'conc_ug_l': np.float32,
'conc_mg_m3': np.float32,
'filtered_vol_l': np.float32,
'vol_after_melt_l': np.float32,
'ice_added_vol_of_water': np.float32,
'ice_total_vol_water': np.float32,
'ice_calculated_vol_ice': np.float32,
'dilution_factor': np.float32,
'quality_control_peaks': 'category',
'nb_core': np.float32,
'nb_core_9cm_dia': np.float32,
'core_area_9cm_dia': np.float32,
'nb_core_14cm_dia': np.float32,
'core_area_14cm_dia': np.float32,
'area_m2': np.float32,
'description_of_sample': 'category',
'n': np.float32,
'observations': 'category',
'method': 'category',
'pi': 'category'
},
parse_dates=['date','date_time']
)
.rename(
columns=dict(
depth_m='depth',
conc_ug_l='chlfluo',
conc_mg_m3='chlHPLC'
)
)
)
df['doy']=df['date'].apply(lambda x:x.timetuple().tm_yday)
return df
chlHPLC_df=load_pigment()
Observations Qikiqtarjuaq 2016: chlorophyll a: select¶
We select the observations of chlorophyll a concentrations at the ice camp 2016 measured by HPLC.
chlHPLC_df=chlHPLC_df[(chlHPLC_df.mission=='ice_camp_2016') \
& (chlHPLC_df.sample_type=='water') \
& (chlHPLC_df.pigment=='Total Chlorophyll a') \
& (chlHPLC_df.method=='HPLC')]
remove negative chlorophyll a
chlHPLC_df['chlHPLC'][ chlHPLC_df['chlHPLC']<0 ]
Series([], Name: chlHPLC, dtype: float32)
find duplicates
g=chlHPLC_df.groupby('doy')['depth'].value_counts()
g.where(g>1).dropna()
doy depth
146 1.5 2.0
10.0 2.0
Name: count, dtype: float64
inspect duplicates for doy=146
chlHPLC_oneday_df=chlHPLC_df[['depth','chlHPLC','doy']][(chlHPLC_df.doy==146) | (chlHPLC_df.doy==148)]
chlHPLC_oneday_df
| depth | chlHPLC | doy | |
|---|---|---|---|
| 20357 | 0.5 | 0.1056 | 146 |
| 20358 | 1.5 | 0.0899 | 146 |
| 20359 | 5.0 | 0.0906 | 146 |
| 20360 | 10.0 | 0.0905 | 146 |
| 20361 | 20.0 | 0.1011 | 146 |
| 20362 | 40.0 | 0.0633 | 146 |
| 20363 | 0.5 | 0.0971 | 148 |
| 20364 | 1.5 | 0.1483 | 146 |
| 20365 | 5.0 | 0.1086 | 148 |
| 20366 | 10.0 | 0.0890 | 146 |
| 20367 | 20.0 | 0.0735 | 148 |
| 20368 | 40.0 | 0.0399 | 148 |
By human judgement, I replace the doy of rows 20364 and 20366 with 148.
chlHPLC_df.loc[[20364,20366],'doy']=148
chlHPLC_oneday_df=chlHPLC_df[['depth','chlHPLC','doy']][(chlHPLC_df.doy==146) | (chlHPLC_df.doy==148)]
chlHPLC_oneday_df
| depth | chlHPLC | doy | |
|---|---|---|---|
| 20357 | 0.5 | 0.1056 | 146 |
| 20358 | 1.5 | 0.0899 | 146 |
| 20359 | 5.0 | 0.0906 | 146 |
| 20360 | 10.0 | 0.0905 | 146 |
| 20361 | 20.0 | 0.1011 | 146 |
| 20362 | 40.0 | 0.0633 | 146 |
| 20363 | 0.5 | 0.0971 | 148 |
| 20364 | 1.5 | 0.1483 | 148 |
| 20365 | 5.0 | 0.1086 | 148 |
| 20366 | 10.0 | 0.0890 | 148 |
| 20367 | 20.0 | 0.0735 | 148 |
| 20368 | 40.0 | 0.0399 | 148 |
Identify days of year for which the number of observations is less than or equal to 2. There are not measurements at enough depths to be able to calculate a meaningful vertical integration.
g=chlHPLC_df.groupby('doy')['doy'].value_counts()
g.where(g<=2).dropna()
doy 120 2.0 Name: count, dtype: float64
chlHPLC_oneday_df=chlHPLC_df[['depth','chlHPLC','doy']][(chlHPLC_df.doy==120)]
chlHPLC_oneday_df
| depth | chlHPLC | doy | |
|---|---|---|---|
| 20299 | 1.5 | 0.0355 | 120 |
| 20300 | 5.0 | 0.0253 | 120 |
By human judgement, I drop these observations.
chlHPLC_df=chlHPLC_df.drop(chlHPLC_df[(chlHPLC_df.doy==120)].index)
Observations Qikiqtarjuaq 2016: chlorophyll a: vertically integrated chlorophyll a (0-100m)¶
in $\mathrm{ mg\ Chla\ m^{-2} }$
#doy: 123
depths_case1 =np.array([ 1.5, 5, 10, 20, 40])
weights_case1=np.array([ 3.25,4.25, 7.5, 15, 70])
#doys: 125,127
depths_case2 =np.array([0.5, 1.5, 5, 20, 40])
weights_case2=np.array([1, 2.25,9.25, 17.5, 70])
#doys:130,132,134,137,139,144,146,148,151,153,155,158,160,162,165,167,169,172,176
depths_case3 =np.array([0.5,1.5, 5, 10, 20, 40])
weights_case3=np.array([1, 2.25, 4.25, 7.5, 15, 70])
#doy: 141
depths_case4 =np.array([0.5,1.5, 5, 10, 40])
weights_case4=np.array([1, 2.25, 4.25,17.5, 75])
#doy: 174
depths_case5 =np.array([0.5,1.5, 10, 20, 40])
weights_case5=np.array([1, 4.75, 9.25,15, 70])
#doys: 179,181,183,186,188,190
depths_case6 =np.array([0.5,1.5, 10, 20,30, 60])
weights_case6=np.array([1, 4.75, 9.25,10,20, 55])
#doys: 193,195,200
depths_case7 =np.array([ 1.5, 10, 20,30, 60, 75])
weights_case7=np.array([ 5.75, 9.25,10,20, 22.5,32.5])
#doys: 197,204
depths_case8 =np.array([ 1.5, 10, 20,30, 45,60])
weights_case8=np.array([ 5.75, 9.25,10,12.5,15,47.5])
nT=365
array1d_iT_obsvintchlHPLC=np.full(nT,np.NaN)
for iT in np.arange(0,nT):
obsvintchlHPLC=np.NaN
chlHPLC=[]
doy=iT+1
chlHPLC_oneday_df=chlHPLC_df[['depth','chlHPLC']][(chlHPLC_df.doy==doy)]
depths=chlHPLC_oneday_df['depth']
if set(depths_case1).issubset(set(depths)):
for depth in depths_case1:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case1)
elif set(depths_case2).issubset(set(depths)):
for depth in depths_case2:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case2)
elif set(depths_case3).issubset(set(depths)):
for depth in depths_case3:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case3)
elif set(depths_case4).issubset(set(depths)):
for depth in depths_case4:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case4)
elif set(depths_case5).issubset(set(depths)):
for depth in depths_case5:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case5)
elif set(depths_case6).issubset(set(depths)):
for depth in depths_case6:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case6)
elif set(depths_case7).issubset(set(depths)):
for depth in depths_case7:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case7)
elif set(depths_case8).issubset(set(depths)):
for depth in depths_case8:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsvintchlHPLC=np.dot(chlHPLC,weights_case8)
array1d_iT_obsvintchlHPLC[iT]=obsvintchlHPLC
Observations Qikiqtarjuaq 2016: chlorophyll $a$: mean chlorophyll a (0-40m)¶
in $mg\ Chla\ m^{-3}$
#doy: 123
depths_case1 =np.array([ 1.5, 5, 10, 20, 40])
weights_case1=np.array([ 3.25, 4.25, 7.5, 15, 10])
#doys: 125,127
depths_case2 =np.array([ 0.5, 1.5, 5, 20, 40])
weights_case2=np.array([ 1, 2.25, 9.25, 17.5, 10])
#doys:130,132,134,137,139,144,146,148,151,153,155,158,160,162,165,167,169,172,176
depths_case3 =np.array([ 0.5, 1.5, 5, 10, 20, 40])
weights_case3=np.array([ 1, 2.25, 4.25, 7.5, 15, 10])
#doy: 141
depths_case4 =np.array([ 0.5, 1.5, 5, 10, 40])
weights_case4=np.array([ 1, 2.25, 4.25, 17.5, 15])
#doy: 174
depths_case5 =np.array([ 0.5, 1.5, 10, 20, 40])
weights_case5=np.array([ 1, 4.75, 9.25, 15, 10])
#doys: 179,181,183,186,188,190
depths_case6 =np.array([ 0.5, 1.5, 10, 20, 30, 60])
weights_case6=np.array([ 1, 4.75, 9.25, 10, 15, 0])
#doys: 193,195,200
depths_case7 =np.array([ 1.5, 10, 20, 30, 60, 75])
weights_case7=np.array([ 5.75, 9.25, 10, 15, 0, 0])
#doys: 197,204
depths_case8 =np.array([ 1.5, 10, 20, 30, 45, 60])
weights_case8=np.array([ 5.75, 9.25, 10, 12.5, 2.5, 0])
nT=365
array1d_iT_obsmeanchlHPLC=np.full(nT,np.NaN)
for iT in np.arange(0,nT):
obsmeanchlHPLC=np.NaN
chlHPLC=[]
doy=iT+1
chlHPLC_oneday_df=chlHPLC_df[['depth','chlHPLC']][(chlHPLC_df.doy==doy)]
depths=chlHPLC_oneday_df['depth']
if set(depths_case1).issubset(set(depths)):
for depth in depths_case1:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case1)
elif set(depths_case2).issubset(set(depths)):
for depth in depths_case2:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case2)
elif set(depths_case3).issubset(set(depths)):
for depth in depths_case3:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case3)
elif set(depths_case4).issubset(set(depths)):
for depth in depths_case4:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case4)
elif set(depths_case5).issubset(set(depths)):
for depth in depths_case5:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case5)
elif set(depths_case6).issubset(set(depths)):
for depth in depths_case6:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case6)
elif set(depths_case7).issubset(set(depths)):
for depth in depths_case7:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case7)
elif set(depths_case8).issubset(set(depths)):
for depth in depths_case8:
chlHPLCtempo \
=chlHPLC_oneday_df['chlHPLC'][(chlHPLC_oneday_df.depth==depth)]
assert chlHPLCtempo.size==1
chlHPLCtempo=chlHPLCtempo.iat[0]
chlHPLC.append(chlHPLCtempo)
obsmeanchlHPLC=np.average(np.array(chlHPLC),weights=weights_case8)
array1d_iT_obsmeanchlHPLC[iT]=obsmeanchlHPLC
Observations Qikiqtarjuaq 2016: chlorophyll $a$: timing¶
Observations Qikiqtarjuaq 2016: chlorophyll $a$: timing: date of bloom initiation¶
The bloom initiation was identified when the accumulation rates ($d^{-1}$) were positive for one week (Boss and Behrenfeld, 2010):
accumulation_rate=$\frac{chl_{iT}-chl_{iT-1}}{chl_{iT-1}}$
$In\ situ$ data was smoothed with LOWESS (locally weighted scatterplot smoothing) and interpolated before calculating the accumulation rates.
For $in\ situ$ data, by applying the definition of Boss & Behrenfeld (2010) of 7 days with a positive accumulation rate ($d^{-1}$), the date of bloom initiation would have been May 9. This date was questionable because light was low on May 9, there was no continuous decrease of nutrients and there was no continuous increase of Chl $a$. Also, changing the number of days with a positive accumulation rate from 7 days to 8 days changed the date of bloom initiation from May 9 to May 27 indicating a high sensitivity on the number of days with a positive accumulation rate. This was an indication that human judgment was needed to identify the date of the bloom initiation. We supposed that the small peak observable on May 9 may not be due to local phytoplankton but instead to the flushing of sea ice algae or to the horizontal advection of phytoplankton. The algorithm was applied after May 9 and found a date of bloom initiation on May 27. Note that Oziel et al. (2019) in $Elem.\ Sci\ Anth.$ also found a date of bloom initiation at the end of May from the same observations.
x=first_year
y=array1d_iT_obsvintchlHPLC
df=pd.DataFrame({'Day of year':x,'Observed Chl a vert. integrated (0-100m)':y})
# Create a scatter plot
fig = px.scatter(df, x=df['Day of year'], y=df['Observed Chl a vert. integrated (0-100m)'],
opacity=0.8, color_discrete_sequence=['black'])
# Change chart background color
fig.update_layout(dict(plot_bgcolor = 'white'))
# Update axes lines
fig.update_xaxes(showgrid=True, gridwidth=1, gridcolor='lightgrey',
zeroline=True, zerolinewidth=1, zerolinecolor='lightgrey',
showline=True, linewidth=1, linecolor='black')
fig.update_yaxes(showgrid=True, gridwidth=1, gridcolor='lightgrey',
zeroline=True, zerolinewidth=1, zerolinecolor='lightgrey',
showline=True, linewidth=1, linecolor='black',
range=[0,251])
# Set figure title
fig.update_layout(title=dict(text="Observed Chl a vert. integrated (0-100m)",
font=dict(color='black')))
# Update marker size
fig.update_traces(marker=dict(size=3))
fig.show()
Figure. Scatterplot of day of year and observed Chl a vert. integrated (0-100m).
y_hat=statsmodels.nonparametric.smoothers_lowess.lowess(y,x,frac=1/4)
# Create a scatter plot
fig = px.scatter(df, x=df['Day of year'], y=df['Observed Chl a vert. integrated (0-100m)'],
opacity=0.8, color_discrete_sequence=['black'])
# Add the prediction line
fig.add_traces(go.Scatter(x=y_hat[:,0], y=y_hat[:,1], name='LOWESS, frac=1/4', line=dict(color='red')))
# Change chart background color
fig.update_layout(dict(plot_bgcolor = 'white'))
# Update axes lines
fig.update_xaxes(showgrid=True, gridwidth=1, gridcolor='lightgrey',
zeroline=True, zerolinewidth=1, zerolinecolor='lightgrey',
showline=True, linewidth=1, linecolor='black')
fig.update_yaxes(showgrid=True, gridwidth=1, gridcolor='lightgrey',
zeroline=True, zerolinewidth=1, zerolinecolor='lightgrey',
showline=True, linewidth=1, linecolor='black',
range=[0,251])
# Set figure title
fig.update_layout(title=dict(text="Observed Chl a vert. integrated (0-100m)",
font=dict(color='black')))
# Update marker size
fig.update_traces(marker=dict(size=3))
fig.show()
# ------- Define interploation functions -------
# Linear - draws a line between the two nearest points and calculates y value based on the slope of that line
f_linear = scipy.interpolate.interp1d(y_hat[:,0],
y=y_hat[:,1],
bounds_error=False)
# Create a new set of points with x values
indices_not_nan=np.argwhere(~np.isnan(y)).flatten()
start=indices_not_nan[0]
end =indices_not_nan[-1]
xnew=np.empty(end-start+1)
for iT in range(start,end+1):
xnew[iT-start]=iT
# Find y values based on the interpolation method
array1d_iT_obsvintchlHPLClowess=f_linear(xnew)
start_edge=np.empty(start)
start_edge[:]=np.NaN
end_edge=np.empty(365-1-end)
end_edge[:]=np.NaN
array1d_iT_obsvintchlHPLClowess=np.concatenate((start_edge,
array1d_iT_obsvintchlHPLClowess,
end_edge))
x=first_year
y=array1d_iT_obsvintchlHPLClowess
df=pd.DataFrame({'Day of year':x,'Smoothed and interpolated observed Chl a vert. integrated (0-100m)':y})
indices_not_nan=np.argwhere(~np.isnan(array1d_iT_obsvintchlHPLClowess)).flatten()
start=indices_not_nan[0]
end =indices_not_nan[-1]
# no NaN between first and last smoothed observed value
assert np.array_equal(indices_not_nan,np.arange(start,end+1)), (
'There is at least one NaN in array1d_iT_obsvintchlHPLCsmooth before '
'its first and last valid values.')
array1d_iT_r2016=np.empty(365)
array1d_iT_r2016[:]=np.NaN
for iT in range(start+1,end):
chliT=array1d_iT_obsvintchlHPLClowess[iT]
chljT=array1d_iT_obsvintchlHPLClowess[iT-1]
acc=(chliT-chljT)/chljT
array1d_iT_r2016[iT]=acc
x=first_year
y=array1d_iT_r2016
df=pd.DataFrame({'Day of year':x,'Accumulation rates':y})
found=False
iT=136 # starting after the small peak of 129 (May 9) + 7
while not found and iT<=364:
array1d_iT_accweekbefore=array1d_iT_r2016[np.arange(iT-6,iT+1)]
if(np.all(~np.isnan(array1d_iT_accweekbefore)) \
and \
np.all(array1d_iT_accweekbefore>0)):
found=True
else:
iT=iT+1
tI=iT
year=2016 # any leap year
# transform day of year into date
# https://stackoverflow.com/questions/2427555/python-question-year-and-day-of-year-to-date
# see answer of Ignacio Vazquez-Abrams
tIs=datetime.datetime(year,1,1) + datetime.timedelta(tI)
Day of year
tI=tI+1
tI
148
tIs.strftime("%B %d")
'May 27'
cI=array1d_iT_obsvintchlHPLClowess[tI-1]
(in $mg\ Chl\ m^{-2}$)
cI
6.4125088304456135
Observations Qikiqtarjuaq 2016: chlorophyll $a$: timing: date of bloom peak¶
The date of the bloom peak was the date of the maximum (vertically integrated) chlorophyll $a$.
tP=np.nanargmax(array1d_iT_obsvintchlHPLC)
tP=tP.item()
year=2016
# transform day of year into date
# https://stackoverflow.com/questions/2427555/python-question-year-and-day-of-year-to-date
# see answer of Ignacio Vazquez-Abrams
tPs=datetime.datetime(year,1,1) + datetime.timedelta(tP)
(in $mg\ Chl\ m^{-2}$)
cP=np.nanmax(array1d_iT_obsvintchlHPLC)
cP
204.44257855415344
Day of year
tP=tP+1
tP
197
tPs.strftime("%B %d")
'July 15'
NEMO-LIM3¶
NEMO-LIM3: simulated sea ice freeze-up¶
infile_sicFile='data/DataS4_output_nemo_lim3/siarea.nemo.2016.365.32bits.bin'
array1d_iT_forsic=np.fromfile(infile_sicFile,dtype='>f4')
The simulated sea ice freeze-up was on 24 November, doy (0-based)=328.
sim_freezeup\
=next(i for i,v in reversed(list(enumerate(array1d_iT_forsic))) if v < 0.50)
sim_freezeup
328
Plot: Observations Qikiqtarjuaq 2016: chlorophyll a¶
Plot in situ phytoplankton accumulation rates ($\mathrm{ d^{-1} }$).
def make_plots(axs):
box=dict(pad=5,alpha=0)
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
xlims=(first_year[0],first_year[-1])
xmin=xlims[0]
xmax=xlims[1]
# --- SNOW AND ICE THICKNESS
ax=axs[0]
ymin=-1.7
ymax=1.0
xsnow=(snow_thickness_df['doy']-1).to_numpy()
ysnow=snow_thickness_df['sample_thickness_m_average'].to_numpy()
h1=ax.plot(xsnow,ysnow,
'o',color='black',markersize=6,label='obs')
xsnow_to_fill=np.arange(xsnow[0],xsnow[-1]+1)
ysnow_to_fill=np.interp(xsnow_to_fill,xsnow,ysnow)
h2=ax.fill_between(xsnow_to_fill,ysnow_to_fill,
linestyle='-',color='aliceblue')
ax.text(95,0.3,'snow',
horizontalalignment = 'left',verticalalignment = 'center')
xice=(ice_thickness_df['doy']-1).to_numpy()
yice=-ice_thickness_df['sample_thickness_m_average'].to_numpy()
h3=ax.plot(xice,yice,'o',color='black',markersize=6)
xice_to_fill=np.arange(xice[0],xice[-1]+1)
yice_to_fill=np.interp(xice_to_fill,xice,yice)
h4=ax.fill_between(xice_to_fill,yice_to_fill,
linestyle='-',color='dodgerblue')
ax.text(95,-0.40,'ice',
horizontalalignment = 'left',verticalalignment = 'center')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylim(ymin,ymax)
ax.yaxis.set_ticks_position('both')
ax.set_ylabel('Thickness\n(m)',bbox=box)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.03,0.75,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=25)
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=364-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
ax.axvline(tI-1,color='green',linestyle=':',label='bloom initiation')
ax.axvline(obs_nosnow,color='black',linewidth=1,linestyle='--')
# --- ACCUMULATION RATE
ax=axs[1]
ymin=-0.5
ymax= 0.5
h1=ax.plot(first_year,array1d_iT_r2016,'-',color='black',lw=3)
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylim(ymin,ymax)
ax.yaxis.set_ticks_position('both')
ax.set_ylabel('Accumulation rate\n($\mathrm{ d^{-1} }$)')
ax.grid()
ax.plot(xlims,[0,0],'k--',color='gray')
plt.text(0.03,0.93,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=25)
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=364-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
ax.axvline(obs_nosnow,color='black',linewidth=1,linestyle='--')
ax.axvline(tI-1,color='green',linestyle=':',label='bloom initiation')
ax.legend(loc='upper right')
# --- ADDITIONAL X-AXIS
ax=axs[2]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0].set_position( [0.12,0.75 ,0.75,0.179])
axs[1].set_position( [0.12,0.107 ,0.75,0.607])
axs[2].set_position( [0.12,0.061 ,0.75,0.032])
with plt.style.context('mplstyles/r.mplstyle'):
# Plot
fig,axs=plt.subplots(3,1,figsize=(18.1,7.858))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_r.'+ext
plt.savefig(outfile,dpi=dpi)
/var/folders/vq/5y09539n705gpb79d6ysx09r0000gr/T/ipykernel_78452/1215027379.py:70: UserWarning: color is redundantly defined by the 'color' keyword argument and the fmt string "k--" (-> color='k'). The keyword argument will take precedence.
Observations from literature¶
obsphotofile='data/DataS7_literature.csv'
obsphoto_df=pd.read_csv(obsphotofile)
Maximum photosynthesis rate at saturation irradiance in nutrient replete conditions ($P^C_{max} \gamma^T$, $\mathrm{ d^{-1} }$).
obsPCmax_df=obsphoto_df[obsphoto_df.parameter=='PCmax']
Lineare initial slope of the Chl a-specific photosynthesis versus irradiance curve ($\alpha^{chl}$, $\mathrm{ g\ C\ (g\ Chl\ a)^{-1}\ h^{-1}\ (\mu mol\ photons\ m^{-2}\ s^{-1})^{-1} }$).
obsalphachl_df=obsphoto_df[obsphoto_df.parameter=='alphastar']
Light saturation parameter ($E_k$, $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$).
obsEk_df=obsphoto_df[obsphoto_df.parameter=='Ek']
Model MITgcm¶
The file car.0000000000.t001.nc was the output of the reference (EXP-0) simulation.
carfile='data/DataS6_output_mitgcm/exp0/car.0000000000.t001.nc'
Model MITgcm: PAR¶
($\mathrm{ mol\ photons\ m^{-2}\ d^{-1} }$)
stod=86400 # 3600*24 s per d
array2d_idepth_iT_PARfullumolEonm2s\
=netcdf_tools.read_netcdf(carfile,'PAR').squeeze().transpose()
if np.all(array2d_idepth_iT_PARfullumolEonm2s==0):
#define GUD_PARUICE
array2d_idepth_iT_PAR_icefullumolEonm2s\
=netcdf_tools.read_netcdf(carfile, 'PAR_ice').squeeze().transpose()
array2d_idepth_iT_PAR_icefull\
=array2d_idepth_iT_PAR_icefullumolEonm2s*1e-6*stod
array2d_idepth_iT_PAR_owfullumolEonm2s\
=netcdf_tools.read_netcdf(carfile, 'PAR_ow').squeeze().transpose()
array2d_idepth_iT_PAR_owfull\
=array2d_idepth_iT_PAR_owfullumolEonm2s*1e-6*stod
(ndepths,ndays)=array2d_idepth_iT_PAR_icefull.shape
array2d_idepth_iT_icefull\
=np.tile(array1d_iT_sic,(ndepths,ndays//365))
array2d_idepth_iT_meanPARfull\
=array2d_idepth_iT_PAR_owfull*(1-array2d_idepth_iT_icefull)\
+array2d_idepth_iT_PAR_icefull*array2d_idepth_iT_icefull
else:
#undef GUD_PARUICE
array2d_idepth_iT_meanPARfull\
=array2d_idepth_iT_PARfullumolEonm2s*1e-6*stod
array2d_idepth_iT_meanPAR=array2d_idepth_iT_meanPARfull[:,array1d_iT1y_iT]
($\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$)
array2d_idepth_iT_meanPARumolEonm2s=array2d_idepth_iT_meanPAR/1e-6/stod
Model MITgcm: nutrients¶
array2d_idepth_iT_modno3full=netcdf_tools.read_netcdf(carfile, 'TRAC04').squeeze().transpose()
array2d_idepth_iT_modno3=array2d_idepth_iT_modno3full[:,array1d_iT1y_iT]
array2d_idepth_iT_modsioh4full=netcdf_tools.read_netcdf(carfile, 'TRAC06').squeeze().transpose()
array2d_idepth_iT_modsioh4=array2d_idepth_iT_modsioh4full[:,array1d_iT1y_iT]
array2d_idepth_iT_modpo4full=netcdf_tools.read_netcdf(carfile, 'TRAC05').squeeze().transpose()
array2d_idepth_iT_modpo4=array2d_idepth_iT_modpo4full[:,array1d_iT1y_iT]
Model MITgcm: Chl $a$¶
indir_exp0='data/DataS6_output_mitgcm/exp0'
chlfile=os.path.join(indir_exp0,'chl.0000000000.t001.nc')
# mg Chl m^-3
array2d_idepth_iT_modchlfull\
=read_mitgcm.get_array2d_idepth_iT_chlfull(chlfile)
array2d_idepth_iT_modchl\
=array2d_idepth_iT_modchlfull[:,array1d_iT1y_iT]
Model MITgcm: Chl $a$: mean Chl $a$ (0-40m)¶
Mean Chl $a$ (0-40m) in $\mathrm{ mg\ Chl\ m^{-3} }$.
array1d_iT_modmeanchl=vstats_pos.mean(
array2d_idepth_iT_tracer=array2d_idepth_iT_modchl,
array1d_idepth_delR=drF,
depth_end=40)
Model MITgcm: $\mathrm{ C_{phyto} }$¶
indir_exp0='data/DataS6_output_mitgcm/exp0'
cphytofile=os.path.join(indir_exp0,'car.0000000000.t001.nc')
# mg C m^-3
array2d_idepth_iT_modcphytofull\
=read_mitgcm.get_array2d_idepth_iT_phytobiofull(cphytofile)
array2d_idepth_iT_modcphyto\
=array2d_idepth_iT_modcphytofull[:,array1d_iT1y_iT]
Model MITgcm: $\mathrm{ C_{phyto} }$: mean $\mathrm{ C_{phyto} }$ (0-40m)¶
Mean phytoplankton biomass (0-40m) in $\mathrm{ mg\ C\ m^{-3} }$.
array1d_iT_modmeancphyto=vstats_pos.mean(
array2d_idepth_iT_tracer=array2d_idepth_iT_modcphyto,
array1d_idepth_delR=drF,
depth_end=40)
Model: photosynthetic parameters¶
Temperature
tavefile='data/DataS6_output_mitgcm/exp0/tave.0000000000.t001.nc'
array2d_idepth_iT_Ttavefull=netcdf_tools.read_netcdf(tavefile,'Ttave').squeeze().transpose()
array2d_idepth_iT_Ttavefull[-1,:]=np.nan
array2d_idepth_iT_Ttave=array2d_idepth_iT_Ttavefull[:,array1d_iT1y_iT]
minT=np.nanmin(array2d_idepth_iT_Ttave)
maxT=np.nanmax(array2d_idepth_iT_Ttave)
IPython.display.Markdown(
'''The interval of temperature ($^{\circ} C$) is
$[%f,%f].$'''
%(minT,maxT)
)
The interval of temperature ($^{\circ} C$) is $[-1.755564,2.919083].$
$\gamma^T$
The modification of growth rate by temperature (unitless, between 0 and 1) is
$\gamma^T=\tau_T e^{A_T(\frac{1}{T+273.15}-\frac{1}{T_N})}$
where
$\tau_T$ is the normalization factor for temperature function $=0.8$ (unitless)
$A_T$ is a coefficient $=-4000$ K
$T$ is the ambient temperature (K)
$T_N$ is the reference temperature $=20\mathrm{ ^{\circ}C }+273.15\mathrm{ ^{\circ}C }=293.15$ K
tau_T=0.8 # unitless
A_T=-4000 # K
T_N=20+273.15 # K
mingammaT=tau_T*math.exp(A_T*(1/(minT+273.15)-1/T_N))
maxgammaT=tau_T*math.exp(A_T*(1/(maxT+273.15)-1/T_N))
IPython.display.Markdown(
'''The interval of $\gamma^T$ (unitless)
in the simulation is
$[%f,%f].$'''
%(mingammaT,maxgammaT)
)
The interval of $\gamma^T$ (unitless) in the simulation is $[0.267952,0.343909].$
Maximum photosynthesis rate at saturation irradiance in nutrient replete conditions ($P^C_{max} \gamma^T$, $\mathrm{ d^{-1} }$)
minpcmaxgammaT=pcmax[15-1]*mingammaT
maxpcmaxgammaT=pcmax[15-1]*maxgammaT
IPython.display.Markdown(
'''The interval of $P^C_{max}\ \gamma^T(T)$ for the diatom 6.6 $\mu$m in $d^{-1}$ is
$[%f,%f].$'''
%(minpcmaxgammaT,maxpcmaxgammaT)
)
The interval of $P^C_{max}\ \gamma^T(T)$ for the diatom 6.6 $\mu$m in $d^{-1}$ is $[0.699651,0.897982].$
Lineare initial slope of the Chl a-specific photosynthesis versus irradiance curve ($\alpha^{chl}$, $\mathrm{ g\ C\ (g\ Chl\ a)^{-1}\ h^{-1}\ (\mu mol\ photons\ m^{-2}\ s^{-1})^{-1} }$)
alphachl=info.f['alphachl'] # mmol C (mg Chl a)^-1 (umol photons m^-2)^-1
alphachl=alphachl[15-1]
For your information, $\alpha^{chl}$ is the same for all numerical phytoplankton types.
We convert the units to compare with literature.
modalphachl_modunits=alphachl # mmol C (mg Chl a)^-1 (umol photons m^-2)^-1
mmolC2molC=1E-3 # mol C (mmol C )^-1
molarmassC=12.011 # g C (mol C)^-1
mgChla2gChla=1E3 # mg Chl a (g Chl a)^-1
s2h=3600 # sh^-1
modalphachl=modalphachl_modunits*mmolC2molC*molarmassC*mgChla2gChla*s2h
IPython.display.Markdown(
r'''$\alpha^{chl} = %.3G
\ g\ C\ (g\ Chl\ a)^{-1}\ h^{-1}\ (\mu mol\ photons\ m^{-2}\ s^{-1})^{-1}$.'''
%(modalphachl)
)
$\alpha^{chl} = 0.0346 \ g\ C\ (g\ Chl\ a)^{-1}\ h^{-1}\ (\mu mol\ photons\ m^{-2}\ s^{-1})^{-1}$.
Light saturation parameter when there is light or no light ($E_k$, $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$)
Maximum photosynthesis rate at saturation irradiance in nutrient replete conditions for numerical diatom 6.6 $\mu$m in $\mathrm{ d^{-1} }$
$P^C_m(j=15) = P^C_{max}(j=15) * \gamma^T(T)$
array2d_idepth_iT_modPCmrep15=pcmax[15-1] \
* tau_T*np.exp(A_T*(1/(array2d_idepth_iT_Ttave+273.15)-1/T_N))
Chlorophyll $a$ for numerical diatom 6.6 $\mu$m ($\mathrm{ mg\ Chl\ a\ m^{-3} }$)
chlfile='data/DataS6_output_mitgcm/exp0/chl.0000000000.t001.nc'
array2d_idepth_iT_diatom7umchlfull\
=netcdf_tools.read_netcdf(chlfile, 'TRAC84').squeeze().transpose()
array2d_idepth_iT_diatom7umchlfull[-1,:]=np.nan
array2d_idepth_iT_diatom7umchl\
=array2d_idepth_iT_diatom7umchlfull[:,array1d_iT1y_iT]
Biomass for numerical diatom 6.6 $\mu$m ($\mathrm{ mmol\ C\ m^{-3} }$)
array2d_idepth_iT_diatom7umbiofull\
=netcdf_tools.read_netcdf(carfile, 'TRAC35').squeeze().transpose()
array2d_idepth_iT_diatom7umbiofull[-1,:]=np.nan
array2d_idepth_iT_diatom7umbio\
=array2d_idepth_iT_diatom7umbiofull[:,array1d_iT1y_iT]
Chl $a$: C ratio for numerical diatom 6.6 $\mu$m ($\mathrm{ mg\ Chl\ {a}\ (mmol\ C)^{-1} }$)
array2d_idepth_iT_modtheta15\
=array2d_idepth_iT_diatom7umchl/array2d_idepth_iT_diatom7umbio
Light saturation parameter for numerical diatom 6.6 $\mu$m ($E_k$, $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$)
$E_k= \frac{P^C_m}{\alpha^{chl} \cdot \theta \cdot \frac{86400\ s}{d}}$
where
$E_k$: light saturation parameter ($\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$).
$P^C_m$: carbon-specific maximum photosynthesis rate ($\mathrm{ d^{-1} }$).
$\alpha^{chl}$: linear initial slope of the Chl a-specific photosynthesis versus irradiance curve ($\mathrm{ mmol\ C\ (mg\ Chl\ a)^{-1}\ (\mu mol\ photons\ m^{-2})^{-1} }$).
$\theta$: Chl a : C ratio ($\mathrm{ mg\ Chl\ a\ (mmol\ C)^{-1} }$).
s2d=86400 # s d^-1
array2d_idepth_iT_modEk15=\
array2d_idepth_iT_modPCmrep15 \
/ (modalphachl_modunits * array2d_idepth_iT_modtheta15 * s2d)
Light saturation parameter when there is light ($E_k$, $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$)
array2d_idepth_iT_modEk15light=array2d_idepth_iT_modEk15.copy()
array2d_idepth_iT_modEk15light[array2d_idepth_iT_meanPARumolEonm2s<5]=np.nan
Plot: compare traits from literature data and simulation¶
obsPCmax_df=obsPCmax_df.copy() # to avoid a SettingWithCopyWarning
obsPCmax_df.sort_values(by=['species'],inplace=True)
obsalphachl_df=obsalphachl_df.copy() # to avoid a SettingWithCopyWarning
obsalphachl_df.sort_values(by=['species'],inplace=True)
obsEk_df=obsEk_df.copy() # to avoid a SettingWithCopyWarning
obsEk_df.sort_values(by=['species'],inplace=True)
array1d_i_modEk15light=array2d_idepth_iT_modEk15light.flatten()
pos=~np.isnan(array1d_i_modEk15light)
array1d_i_modEk15light=array1d_i_modEk15light[pos]
def make_plots(axs):
# --- MAXIMUM PHOTOSYNTHESIS RATE AT SATURATION IRRADIANCE
# --- IN NUTRIENT REPLETE CONDITIONS
# --- (P^C_max * gamma^T)
ax=axs[0,0]
xlim=(0,1.5)
h1=ax.scatter(obsPCmax_df.value,obsPCmax_df.species,
color='b',label='obs')
h2=ax.axhline('',
xmin=minpcmaxgammaT/xlim[1],
xmax=maxpcmaxgammaT/xlim[1],
color='red',lw=10)
h3=ax.axhline('foo',
lw=0)
# to have the ytick labels in alphabetical order from top to bottom
ax.invert_yaxis()
# in italics
# FULL GENUS
ax.set_yticklabels(['Chaetoceros furcellatus',
'Chaetoceros neogracile',
'Fragilariopsis cylindrus',
'Thalassiosira gravida 1',
'Thalassiosira gravida 2',
'Thalassiosira nordenskioeldii',
'Numerical diatom 6.6 μm'
],
fontstyle='italic',
fontweight='bold')
ax.set_xlim(xlim)
ax.set_xlabel('$P^C_{max} \gamma^T$\n($\mathrm{ d^{-1} }$)')
ax.grid()
plt.text(0.10,0.88,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.legend([h1,h2],['obs','model'],ncol=2,loc='lower center')
# --- NOTHING
axs[1,0].axis('off')
# --- LINEAR INITIAL SLOPE OF THE CHL A-SPECIFIC PHOTOSYNTHESIS
# --- VERSUS IRRADIANCE CURVE
# --- (ALPHA^CHL)
ax=axs[0,1]
xlim=(0,0.05)
h1=ax.scatter(obsalphachl_df.value,
obsalphachl_df.species,
c='b',label='obs')
h2=ax.scatter(modalphachl,
'',
color='red')
h3=ax.axhline('foo',
lw=0)
# to have the ytick labels in alphabetical order from top to bottom
ax.invert_yaxis()
ax.set_yticklabels([])
ax.set_xlim(xlim)
ax.set_xlabel('''$\\alpha^{chl}$\n\
$\mathrm{ (g\ C\ (g\ Chl)^{-1}\ h^{-1} }$\n\
$\mathrm{ (\mu mol\ photons\ m^{-2}\ s^{-1})^{-1}) }$''')
ax.grid()
plt.text(0.10,0.88,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
# --- NOTHING
axs[1,1].axis('off')
# --- LIGHT SATURATION PARAMETER
# --- (E_k)
ax=axs[0,2]
xlim=(0,250)
h1=ax.scatter(obsEk_df.value,
obsEk_df.species,
c='b',label='obs')
h3=ax.axhline('foo',
lw=0)
h4=ax.axhline('foo2',
lw=0)
# to have the ytick labels in alphabetical order from top to bottom
ax.invert_yaxis()
ax.set_yticklabels([])
ax.set_xlabel('$E_k$\n$\mathrm{ (\mu mol\ photons\ m^{-2}\ s^{-1}) }$')
ax.set_xlim(xlim)
ax.grid()
plt.text(0.10,0.88,' c ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax=axs[1,2]
boxplot=ax.boxplot(array1d_i_modEk15light,
vert=False,
widths=0.75,
boxprops=dict(color='red'),
flierprops=dict(marker='o',
markerfacecolor='red',
markeredgecolor='none',
markersize=1),
medianprops=dict(color='red'),
capprops=dict(color='red'),
whiskerprops=dict(color='red'))
ax.set_xlim(xlim)
ax.set_yticklabels([])
ax.set_xticklabels([])
ax.spines['top'].set_visible(False)
ax.spines['right'].set_visible(False)
ax.spines['bottom'].set_visible(False)
# ax.spines['left'].set_visible(False)
ax.tick_params(length=0)
ax.grid()
# --- POSITION
# FULL GENUS
axs[0,0].set_position( [0.24,0.21 ,0.22,0.77])
axs[0,1].set_position( [0.50,0.21 ,0.22,0.77])
axs[0,2].set_position( [0.76,0.21 ,0.22,0.77])
axs[1,2].set_position( [0.76,0.31 ,0.22,0.07])
return boxplot
with plt.style.context('mplstyles/lacour17.mplstyle'):
# Plot
fig,axs=plt.subplots(2,3,figsize=(16, 8))
boxplot=make_plots(axs)
# --- SAVE
plt.savefig('figures_timing_sup_progress/lacour17.'+ext,dpi=dpi)
/var/folders/vq/5y09539n705gpb79d6ysx09r0000gr/T/ipykernel_78452/3592232768.py:20: UserWarning: FixedFormatter should only be used together with FixedLocator
Where IQR is the interquartile range (Q3-Q1), the upper whisker will extend to last datum less than Q3+whis*IQR (https://matplotlib.org/3.1.1/api/_as_gen/matplotlib.pyplot.boxplot.html).
upper_whisker=[item.get_xdata()[1] for item in boxplot['whiskers']][1]
fliers=array1d_i_modEk15light[array1d_i_modEk15light>upper_whisker]
nb_fliers=fliers.size
nb_total=array1d_i_modEk15light.size
Number of outliers
nb_fliers
357
Number of data points
nb_total
3357
Output from model CGRF¶
input file: 1D_GDPS_qsw_y2016.nc
name: qsw
Variable qsw is the downwelling shortwave radiation just above surface in energy units ($\mathrm{ E_s(z=0^+,SW) }$ in $\mathrm{ W m^{-2} }$). It is at indices (2,2) meaning centre of grid point was used. File 1D_GDPS_qsw_y2016.nc is the output of CMC GDPS reforecasts (CGRF, Smith et al., 2014, https://doi.org/10.1002/qj.2194.) CMC: Canadian Meteorological Centre's. GDPS: global deterministic prediction system. The data was provided by Gregory Smith (Environment Canada). Gaetan Olivier (UBO) preprocessed it to get the data at the Green Edge 2016 ice camp (67.48N, -63.79E).
The goal was to calculate scalar PAR just below open water in photons density flux ($\mathrm{ E_{0,w}(z=0^-,PAR[Q]) }$ in $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$) from $\mathrm{ E_s(z=0^+,SW) }$.
Daily downwelling shortwave radiation just above surface in energy units ($\mathrm{ E_s(z=0^+,SW) }$ in $\mathrm{ W m^{-2} }$).
ndays=365
qswow_file='data/DataS8_output_cgrf/1D_GDPS_qsw_y2016.nc'
array1d_iT_qswowabovehourly=netcdf_tools.read_netcdf(qswow_file,'qsw')
array1d_iT_qswowabovehourly=array1d_iT_qswowabovehourly[:,1,1]
array1d_idoy_qswowabove=np.zeros(ndays)
for idoy in first_year:
qswowabovedaily=0
for iT in np.arange(idoy*24,(idoy+1)*24):
qswowabovedaily=qswowabovedaily+array1d_iT_qswowabovehourly[iT]
qswowabovedaily=qswowabovedaily/24
array1d_idoy_qswowabove[idoy]=qswowabovedaily
Daily downwelling shortwave radiation just below open water in energy units ($\mathrm{ E_s(z=0^-,SW) }$ in $\mathrm{ W m^{-2} }$).
albedo=0.066 # model PISCES, ref.: Aumont et al. (2015)
array1d_idoy_qswow=array1d_idoy_qswowabove*(1-albedo)
Daily scalar PAR just below open water in energy units ($\mathrm{ E_{0,w}(z=0^-,PAR[W]) }$ in $\mathrm{ W m^{-2} }$).
Wm2_to_PAR=0.43 # model PISCES, ref.: Aumont et al. (2015)
array1d_idoy_PARow=array1d_idoy_qswow*Wm2_to_PAR
Daily scalar irradiance just below open water for band lambda in energy units ($\mathrm{ E_{0,w,\lambda}(z=0^-,PAR[W]) }$ in $\mathrm{ W m^{-2} }$).
onecolor=1/3
array1d_idoy_blueWm2ow=array1d_idoy_PARow*onecolor
array1d_idoy_greenWm2ow=array1d_idoy_PARow*onecolor
array1d_idoy_redWm2ow=array1d_idoy_PARow*onecolor
Daily scalar irradiance just below open water for band lambda in photons density flux ($\mathrm{ E_{0,w,\lambda}(z=0^-,PAR[Q]) }$ in $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$).
rmus=1/0.83
planck = 6.6256e-34 # Plancks constant J sec
c = 2.998e8 # speed of light m/sec
hc = 1.0/(planck*c)
oavo = 1.0/6.023e23 # 1/Avogadros number
hcoavo = hc*oavo
wb_center=np.array([450,550,650]) # [blue, green, red]
rlamm=wb_center*1e-9 # lambda in m
WtouEins=1e6*rlamm*hcoavo
array1d_idoy_blueow =array1d_idoy_blueWm2ow *WtouEins[0]
array1d_idoy_greenow=array1d_idoy_greenWm2ow*WtouEins[1]
array1d_idoy_redow =array1d_idoy_redWm2ow *WtouEins[2]
Daily scalar PAR just below open water in photons density flux ($\mathrm{ E_{0,w}(z=0^-,PAR[Q]) }$ in $\mathrm{ \mu mol\ photons\ m^{-2}\ s^{-1} }$).
array1d_idoy_PARowuEin=array1d_idoy_blueow\
+array1d_idoy_greenow\
+array1d_idoy_redow
Plot: PAR¶
def make_plots(axs):
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
# --- DOWNWELLING SHORTWAVE RADIATION JUST ABOVE SURFACE
# --- IN ENERGY UNITS
ax=axs[0]
ax.plot(first_year,array1d_idoy_qswowabove,color='black')
ax.set_ylabel('$\mathrm{ E_s(z=0^+,SW) }$\n'
'($\mathrm{ W\ m^{-2} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,400)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.20,
'JUST ABOVE SURFACE',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
plt.text(0.50,1.05,
'(a) Shortwave radiation',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- DOWNWELLING SHORTWAVE RADIATION JUST BELOW OPEN WATER
# --- IN ENERGY UNITS
ax=axs[1]
ax.plot(first_year,array1d_idoy_qswow,color='black')
ax.set_ylabel('$\mathrm{ E_s(z=0^-,SW) }$\n'
'($\mathrm{ W\ m^{-2} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,400)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.20,
'JUST BELOW OPEN WATER',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
plt.text(0.50,1.05,
'(b) Shortwave radiation',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- DAILY SCALAR PAR JUST BELOW OPEN WATER IN ENERGY UNITS
ax=axs[2]
ax.plot(first_year,array1d_idoy_PARow,color='black')
ax.set_ylabel('$\mathrm{ E_{0,w}(z=0^-,PAR[W]) }$\n'
'($\mathrm{ W\ m^{-2} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,400)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.05,
'(c) PAR',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- DAILY SCALAR IRRADIANCE JUST BELOW OPEN WATER FOR BAND LAMBDA
# --- IN ENERGY UNITS
ax=axs[3]
ax.plot(first_year,array1d_idoy_blueWm2ow,
':',color='blue',label='blue')
ax.plot(first_year,array1d_idoy_greenWm2ow,
'--',color='green',label='green')
ax.plot(first_year,array1d_idoy_redWm2ow,
'-.',color='red',label='red')
ax.legend()
ax.set_ylabel('$\mathrm{ E_{0,w,\lambda}(z=0^-,PAR[W]) }$\n'
'($\mathrm{ W\ m^{-2} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,400)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.05,
'(d) PAR',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- DAILY SCALAR IRRADIANCE JUST BELOW OPEN WATER FOR BAND LAMBDA
# --- IN PHOTONS DENSITY FLUX
ax=axs[4]
ax.plot(first_year,array1d_idoy_blueow,
':',color='blue',label='blue')
ax.plot(first_year,array1d_idoy_greenow,
'--',color='green',label='green')
ax.plot(first_year,array1d_idoy_redow,
'-.',color='red',label='red')
ax.legend()
ax.set_ylabel('$\mathrm{ E_{0,w,\lambda}(z=0^-,PAR[Q]) }$\n'
'($\mathrm{ \mu\ mol\ photons\ m^{-2}\ s^{-1} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,400)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.05,
'(e) PAR',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- DAILY SCALAR PAR JUST BELOW OPEN WATER IN PHOTONS DENSITY FLUX
ax=axs[5]
ax.plot(first_year,array1d_idoy_PARowuEin,color='black')
ax.set_ylabel('$\mathrm{ E_{0,w}(z=0^-,PAR[Q]) }$\n'
'($\mathrm{ \mu\ mol\ photons\ m^{-2}\ s^{-1} }$)')
ax.set_xlim(first_year[0],first_year[-1])
ax.set_ylim(0,800)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.grid()
plt.text(0.50,1.05,
'(f) PAR',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- ADDITIONAL X-AXIS
ax=axs[6]
ax.set_xlim(first_year[0],first_year[-1])
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0].set_position( [0.13,0.81 ,0.85,0.13])
axs[1].set_position( [0.13,0.63 ,0.85,0.13])
axs[2].set_position( [0.13,0.48 ,0.85,0.13])
axs[3].set_position( [0.13,0.33 ,0.85,0.13])
axs[4].set_position( [0.13,0.18 ,0.85,0.13])
axs[5].set_position( [0.13,0.03 ,0.85,0.13])
axs[6].set_position( [0.13,0.014 ,0.85,0.007])
with plt.style.context('mplstyles/par.mplstyle'):
# Plot
fig,axs=plt.subplots(nrows=7,figsize=(16.00, 35.05))
make_plots(axs)
# --- SAVE
plt.savefig('figures_timing_sup_progress/par.'+ext,dpi=dpi)
Plot: Concentrations of nitrate, silicic acid, phosphate, and Chl $a$¶
def make_plots(axs):
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
RF_above151=RF[RF<151]
# --- NITRATE
ax=axs[0]
h=ax.pcolormesh(first_year366,
RF_above151,
array2d_idepth_iT_modno3[0:(RF_above151.size)-1,:],
cmap='viridis',
vmin=0,
vmax=6)
ax.scatter(nutrients_df['doy']-1,
nutrients_df['depth'],
c=nutrients_df['no3'],
cmap='viridis',
vmin=0,
vmax=6,
s=100,
edgecolor='k')
ax.set_xticklabels([])
ax.set_ylabel('Depth (m)')
ax.set_xlim(first_year366[0],first_year[-1])
ax.set_ylim(150,0)
ax.grid()
ax.set_xticks(locs)
ax.set_xticklabels(labels)
cbaxes = fig.add_axes([0.92,0.78, 0.01,0.19])
cbar=plt.colorbar(h,
cax = cbaxes)
cbar.set_label('(mmol N $\mathrm{ m^{-3} }$)')
plt.text(0.50,1.10,'(a) Nitrate',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- SILICIC ACID
ax=axs[1]
h=ax.pcolormesh(first_year366,
RF_above151,
array2d_idepth_iT_modsioh4[0:(RF_above151.size)-1,:],
cmap='viridis',
vmin=0,
vmax=10)
ax.scatter(nutrients_df['doy']-1,
nutrients_df['depth'],
c=nutrients_df['sioh4'],
cmap='viridis',
vmin=0,
vmax=10,
s=100,
edgecolor='k')
ax.set_xticklabels([])
ax.set_ylabel('Depth (m)')
ax.set_xlim(first_year366[0],first_year[-1])
ax.set_ylim(150,0)
ax.grid()
ax.set_xticks(locs)
ax.set_xticklabels(labels)
cbaxes = fig.add_axes([0.92,0.53, 0.01,0.19])
cbar=plt.colorbar(h,
cax = cbaxes)
cbar.set_label('(mmol Si$\mathrm{ {(OH)}_4 }$ $\mathrm{ m^{-3} }$)')
plt.text(0.50,1.10,'(b) Silicic acid',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- PHOSPHATE
ax=axs[2]
h=ax.pcolormesh(first_year366,
RF_above151,
array2d_idepth_iT_modpo4[0:(RF_above151.size)-1,:],
cmap='viridis',
vmin=0,
vmax=1)
ax.scatter(nutrients_df['doy']-1,
nutrients_df['depth'],
c=nutrients_df['po4'],
cmap='viridis',
vmin=0,
vmax=1,
s=100,
edgecolor='k')
ax.set_xticklabels([])
ax.set_ylabel('Depth (m)')
ax.set_xlim(first_year366[0],first_year[-1])
ax.set_ylim(150,0)
ax.grid()
ax.set_xticks(locs)
ax.set_xticklabels(labels)
cbaxes = fig.add_axes([0.92,0.28, 0.01,0.19])
cbar=plt.colorbar(h,
cax = cbaxes)
cbar.set_label('(mmol $\mathrm{ PO_4 }$ $\mathrm{ m^{-3} }$)')
plt.text(0.50,1.10,'(c) Phosphate',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- CHLOROPHYLL A
ax=axs[3]
h=ax.pcolormesh(first_year366,
RF_above151,
array2d_idepth_iT_modchl[0:(RF_above151.size)-1,:],
cmap='viridis',
norm=mpl.colors.LogNorm(vmin=1E-2, vmax=10))
ax.scatter(chlHPLC_df['doy']-1,
chlHPLC_df['depth'],
c=chlHPLC_df['chlHPLC'],
cmap='viridis',
norm=mpl.colors.LogNorm(vmin=1E-2, vmax=10),
s=100,
edgecolor='k')
ax.set_xticklabels([])
ax.set_ylabel('Depth (m)')
ax.set_xlim(first_year366[0],first_year[-1])
ax.set_ylim(150,0)
ax.grid()
ax.set_xticks(locs)
ax.set_xticklabels(labels)
cbaxes = fig.add_axes([0.92,0.03, 0.01,0.19])
cbar=plt.colorbar(h,
cax = cbaxes)
cbar.set_label('(mg Chl $\mathrm{ m^{-3} }$)')
plt.text(0.50,1.10,'(d) Chlorophyll $a$',
horizontalalignment = 'center', verticalalignment='center',
transform=ax.transAxes)
# --- POSITION
axs[0].set_position( [0.07,0.78 ,0.83,0.19])
axs[1].set_position( [0.07,0.53 ,0.83,0.19])
axs[2].set_position( [0.07,0.28 ,0.83,0.19])
axs[3].set_position( [0.07,0.03 ,0.83,0.19])
with plt.style.context('mplstyles/concentrations_nutrients_Chla.mplstyle'):
# Plot
fig,axs=plt.subplots(4,1,figsize=(16,16))
make_plots(axs)
# --- SAVE
plt.savefig('figures_timing_sup_progress/sup_concentrations_nutrients_Chla.'+ext,dpi=dpi)
EXP-1¶
Differing minimum biomass.
EXP-1: infiles¶
infile_car_EXP1_grp_Xmin_1Eminus3\
='data/DataS6_output_mitgcm/exp1_Cmin1Eminus3/car.0000000000.t001.nc'
infile_car_EXP1_grp_Xmin_1Eminus2\
='data/DataS6_output_mitgcm/exp1_Cmin1Eminus2/car.0000000000.t001.nc'
infile_car_EXP1_grp_Xmin_1Eminus1\
='data/DataS6_output_mitgcm/exp1_Cmin1Eminus1/car.0000000000.t001.nc'
EXP-1: $\mathrm{ \sum C_{phyto} }$¶
$\mathrm{ \sum C_{phyto} }$ is the vertically integrated $\mathrm{ C_{phyto} }$ (0-100m) in $\mathrm{ mg\ C\ m^{-2} }$.
depth_end=100 # m
EXP-1: $\mathrm{ \sum C_{phyto} }$: groups¶
$\mathrm{ \sum C_{phyto} }$ for each phytoplankton group.
EXP-1: $\mathrm{ \sum C_{phyto} }$: groups: grp_Xmin=$\mathrm{ 1*10^{-3}\ mmol\ C\ m^{-3}\ {type}^{-1} }$¶
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biodino_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_bioother_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biopico_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
EXP-1: $\mathrm{ \sum C_{phyto} }$: groups: grp_Xmin=$\mathrm{ 1*10^{-2}\ mmol\ C\ m^{-3}\ {type}^{-1} }$ (default)¶
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biodino_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_bioother_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biopico_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
EXP-1: $\mathrm{ \sum C_{phyto} }$: groups: grp_Xmin=$\mathrm{ 1*10^{-1}\ mmol\ C\ m^{-3}\ {type}^{-1} }$¶
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biodino_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_bioother_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biopico_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
EXP-1: $\mathrm{ \sum C_{phyto} }$: sizes¶
$\mathrm{ \sum C_{phyto} }$ for each phytoplankton size class.
EXP-1: $\mathrm{ \sum C_{phyto} }$: sizes: grp_Xmin=$\mathrm{ 1*10^{-3}\ mmol\ C\ m^{-3}\ {type}^{-1} }$¶
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus3\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP1_grp_Xmin_1Eminus3,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus3\
=array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus3/1000 # mg C -> g C
EXP-1: $\mathrm{ \sum C_{phyto} }$: sizes: grp_Xmin=$\mathrm{ 1*10^{-3}\ mmol\ C\ m^{-2}\ {type}^{-1} }$¶
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus2\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP1_grp_Xmin_1Eminus2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus2\
=array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus2/1000 # mg C -> g C
EXP-1: $\mathrm{ \sum C_{phyto} }$: sizes: grp_Xmin=$\mathrm{ 1*10^{-3}\ mmol\ C\ m^{-1}\ {type}^{-1} }$¶
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus1\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP1_grp_Xmin_1Eminus1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus1\
=array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus1/1000 # mg C -> g C
def make_plots(axs):
box=dict(pad=5,alpha=0)
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
locs=np.array([151,
181, 212, 243, 273, 304, 334])
labels=('Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
xmin=135 # 15 May (0-based)
xmax=364 # 31 December (0-based)
xlims=(xmin,xmax)
ymin=0
ymax=5
ylims=(ymin,ymax)
# --- (A) EXP-1: SUM C_PHYTO: GROUPS: GRP_XMIN=1*10^-3
ax=axs[0,0]
h_biodiatoms_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus3,
'-',color='red',label='diatom')
h_biodino_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus3,
'-',color='magenta',label='dino')
h_bioother_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus3,
'-',color='blue',label='other')
h_biopico_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus3,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-1\n$\mathrm{ C_{min}=10^{-3} }$',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (B) EXP-1: SUM C_PHYTO: GROUPS: GRP_XMIN=1*10^-2 (DEFAULT)
ax=axs[0,1]
h_biodiatoms_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus2,
'-',color='red',label='diatom')
h_biodino_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus2,
'-',color='magenta',label='dino')
h_bioother_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus2,
'-',color='blue',label='other')
h_biopico_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus2,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
plt.text(0.05,0.90,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-1\n$\mathrm{ C_{min}=10^{-2} }$',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (C) EXP-1: SUM C_PHYTO: GROUPS: GRP_XMIN=1*10^-1
ax=axs[0,2]
h_biodiatoms_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP1_grp_Xmin_1Eminus1,
'-',color='red',label='diatom')
h_biodino_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biodino_EXP1_grp_Xmin_1Eminus1,
'-',color='magenta',label='dino')
h_bioother_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_bioother_EXP1_grp_Xmin_1Eminus1,
'-',color='blue',label='other')
h_biopico_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biopico_EXP1_grp_Xmin_1Eminus1,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
plt.text(0.05,0.90,' c ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-1\n$\mathrm{ C_{min}=10^{-1} }$',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (D) EXP-1: SUM C_PHYTO: SIZES: GRP_XMIN=1*10^-3
ax=axs[1,0]
h_biophytolte4_5_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus3,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus3,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus3,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus3,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP1_grp_Xmin_1Eminus3=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus3,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' d ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (E) EXP-1: SUM C_PHYTO: SIZES: GRP_XMIN=1*10^-2
ax=axs[1,1]
h_biophytolte4_5_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus2,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus2,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus2,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus2,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP1_grp_Xmin_1Eminus2=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus2,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' e ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (F) EXP-1: SUM C_PHYTO: SIZES: GRP_XMIN=1*10^-1
ax=axs[1,2]
h_biophytolte4_5_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP1_grp_Xmin_1Eminus1,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP1_grp_Xmin_1Eminus1,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biophyto10_EXP1_grp_Xmin_1Eminus1,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biophyto15_EXP1_grp_Xmin_1Eminus1,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP1_grp_Xmin_1Eminus1=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP1_grp_Xmin_1Eminus1,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' f ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- ADDITIONAL X-AXIS
ax=axs[2,0]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,1]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,2]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0,0].set_position([0.074,0.460, 0.292,0.372])
axs[0,1].set_position([0.388,0.460, 0.292,0.372])
axs[0,2].set_position([0.702,0.460, 0.292,0.372])
axs[1,0].set_position([0.074,0.066, 0.292,0.372])
axs[1,1].set_position([0.388,0.066, 0.292,0.372])
axs[1,2].set_position([0.702,0.066, 0.292,0.372])
axs[2,0].set_position([0.074,0.037, 0.292,0.020])
axs[2,1].set_position([0.388,0.037, 0.292,0.020])
axs[2,2].set_position([0.702,0.037, 0.292,0.020])
with plt.style.context('mplstyles/exp1_groups.mplstyle'):
# Plot
fig,axs=plt.subplots(nrows=3,ncols=3,figsize=(29.445,12.825))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_exp1_groups.'+ext
plt.savefig(outfile,dpi=dpi)
Figure <sup_exp1_groups>. Sensitivity simulations by phytoplankton groups and size classes: EXP-1 prescribed minimum biomass. The minimum biomass ($\mathrm{ C_{min} }$) was prescribed at $10^{-3}$, $10^{-2}$ and $10^{-1}$ $\mathrm{ mmol\ C\ m^{-3} }$ for each of the 26 numerical phytoplankton types. First row: Vertically integrated biomass for each functional group (0-100m). pico: picophytoplankton. other nano: other nanophytoplankton. dino: mixotrophic dinoflagellate. Second row: Vertically integrated biomass (0-100m) for different size classes or intervals of size classes. The size classes were described with their equivalent spherical diameter (ESD). Note that the interval of size classes <= 4.5 μm gathered together four size classes of the model and that the interval of size classes >= 20 μm gathered together eight size classes of the model. The vertical dotted line is the date of the snow melt completion. The grey rectangles represent the time of year before the sea ice break-up and after the sea ice freeze-up.
EXP-2¶
Differing light under sea ice.
EXP-2: infiles¶
infile_car_EXP2_2\
='data/DataS6_output_mitgcm/exp2_2/car.0000000000.t001.nc'
infile_car_EXP2_1\
='data/DataS6_output_mitgcm/exp2_1/car.0000000000.t001.nc'
infile_car_EXP0\
='data/DataS6_output_mitgcm/exp0/car.0000000000.t001.nc'
EXP-2: $\mathrm{ \sum C_{phyto} }$¶
$\mathrm{ \sum C_{phyto} }$ is the vertically integrated $\mathrm{ C_{phyto} }$ (0-100m) in $\mathrm{ mg\ C\ m^{-2} }$.
depth_end=100 # m
EXP-2: $\mathrm{ \sum C_{phyto} }$: groups¶
$\mathrm{ \sum C_{phyto} }$ for each phytoplankton group.
EXP-2: $\mathrm{ \sum C_{phyto} }$: groups: EXP-2.2, opaque snow and opaque ice¶
array1d_iT_biodiatoms_EXP2_2\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP2_2\
=array1d_iT_biodiatoms_EXP2_2/1000 # mg C -> g C
array1d_iT_biodino_EXP2_2\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP2_2\
=array1d_iT_biodino_EXP2_2/1000 # mg C -> g C
array1d_iT_bioother_EXP2_2\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP2_2\
=array1d_iT_bioother_EXP2_2/1000 # mg C -> g C
array1d_iT_biopico_EXP2_2\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP2_2\
=array1d_iT_biopico_EXP2_2/1000 # mg C -> g C
EXP-2: $\mathrm{ \sum C_{phyto} }$: groups: EXP-2.1, opaque snow and translucent ice¶
array1d_iT_biodiatoms_EXP2_1\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP2_1\
=array1d_iT_biodiatoms_EXP2_1/1000 # mg C -> g C
array1d_iT_biodino_EXP2_1\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP2_1\
=array1d_iT_biodino_EXP2_1/1000 # mg C -> g C
array1d_iT_bioother_EXP2_1\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP2_1\
=array1d_iT_bioother_EXP2_1/1000 # mg C -> g C
array1d_iT_biopico_EXP2_1\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP2_1\
=array1d_iT_biopico_EXP2_1/1000 # mg C -> g C
EXP-2: $\mathrm{ \sum C_{phyto} }$: groups: EXP-0, translucent snow and translucent ice (reference)¶
array1d_iT_biodiatoms_EXP0\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP0\
=array1d_iT_biodiatoms_EXP0/1000 # mg C -> g C
array1d_iT_biodino_EXP0\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP0\
=array1d_iT_biodino_EXP0/1000 # mg C -> g C
array1d_iT_bioother_EXP0\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP0\
=array1d_iT_bioother_EXP0/1000 # mg C -> g C
array1d_iT_biopico_EXP0\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP0\
=array1d_iT_biopico_EXP0/1000 # mg C -> g C
EXP-2: $\mathrm{ \sum C_{phyto} }$: sizes¶
$\mathrm{ \sum C_{phyto} }$ for each phytoplankton size class.
EXP-2: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-2.2, opaque snow and opaque ice¶
array1d_iT_biophytolte4_5_EXP2_2\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP2_2\
=array1d_iT_biophytolte4_5_EXP2_2/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP2_2\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP2_2\
=array1d_iT_biophyto6_6_EXP2_2/1000 # mg C -> g C
array1d_iT_biophyto10_EXP2_2\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP2_2\
=array1d_iT_biophyto10_EXP2_2/1000 # mg C -> g C
array1d_iT_biophyto15_EXP2_2\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP2_2\
=array1d_iT_biophyto15_EXP2_2/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP2_2\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP2_2,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP2_2\
=array1d_iT_biophytoge20_EXP2_2/1000 # mg C -> g C
EXP-2: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-2.1, opaque snow and translucent ice¶
array1d_iT_biophytolte4_5_EXP2_1\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP2_1\
=array1d_iT_biophytolte4_5_EXP2_1/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP2_1\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP2_1\
=array1d_iT_biophyto6_6_EXP2_1/1000 # mg C -> g C
array1d_iT_biophyto10_EXP2_1\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP2_1\
=array1d_iT_biophyto10_EXP2_1/1000 # mg C -> g C
array1d_iT_biophyto15_EXP2_1\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP2_1\
=array1d_iT_biophyto15_EXP2_1/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP2_1\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP2_1,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP2_1\
=array1d_iT_biophytoge20_EXP2_1/1000 # mg C -> g C
EXP-2: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-0, translucent snow and translucent ice (reference)¶
array1d_iT_biophytolte4_5_EXP0\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP0\
=array1d_iT_biophytolte4_5_EXP0/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP0\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP0\
=array1d_iT_biophyto6_6_EXP0/1000 # mg C -> g C
array1d_iT_biophyto10_EXP0\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP0\
=array1d_iT_biophyto10_EXP0/1000 # mg C -> g C
array1d_iT_biophyto15_EXP0\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP0\
=array1d_iT_biophyto15_EXP0/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP0\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP0,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP0\
=array1d_iT_biophytoge20_EXP0/1000 # mg C -> g C
def make_plots(axs):
box=dict(pad=5,alpha=0)
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
locs=np.array([151,
181, 212, 243, 273, 304, 334])
labels=('Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
xmin=135 # 15 May (0-based)
xmax=364 # 31 December (0-based)
xlims=(xmin,xmax)
ymin=0
ymax=5
ylims=(ymin,ymax)
# --- (A) EXP-2: SUM C_PHYTO: GROUPS:
# --- EXP-2.2, OPAQUE SNOW AND OPAQUE ICE
ax=axs[0,0]
h_biodiatoms_EXP2_2=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP2_2,
'-',color='red',label='diatom')
h_biodino_EXP2_2=ax.plot(
first_year,
array1d_iT_biodino_EXP2_2,
'-',color='magenta',label='dino')
h_bioother_EXP2_2=ax.plot(
first_year,
array1d_iT_bioother_EXP2_2,
'-',color='blue',label='other')
h_biopico_EXP2_2=ax.plot(
first_year,
array1d_iT_biopico_EXP2_2,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-2.2\nopaque snow\nopaque ice',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (B) EXP-2: SUM C_PHYTO: GROUPS:
# --- EXP-2.1, OPAQUE SNOW AND TRANSLUCENT ICE
ax=axs[0,1]
h_biodiatoms_EXP2_1=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP2_1,
'-',color='red',label='diatom')
h_biodino_EXP2_1=ax.plot(
first_year,
array1d_iT_biodino_EXP2_1,
'-',color='magenta',label='dino')
h_bioother_EXP2_1=ax.plot(
first_year,
array1d_iT_bioother_EXP2_1,
'-',color='blue',label='other')
h_biopico_EXP2_1=ax.plot(
first_year,
array1d_iT_biopico_EXP2_1,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-2.1\nopaque snow\ntranslucent ice',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (C) EXP-2: SUM C_PHYTO: GROUPS:
# --- EXP-0, TRANSLUCENT SNOW AND TRANSLUCENT ICE
ax=axs[0,2]
h_biodiatoms_EXP0=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP0,
'-',color='red',label='diatom')
h_biodino_EXP0=ax.plot(
first_year,
array1d_iT_biodino_EXP0,
'-',color='magenta',label='dino')
h_bioother_EXP0=ax.plot(
first_year,
array1d_iT_bioother_EXP0,
'-',color='blue',label='other')
h_biopico_EXP0=ax.plot(
first_year,
array1d_iT_biopico_EXP0,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' c ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-0\ntranslucent snow\ntranslucent ice',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (D) EXP-2: SUM C_PHYTO: SIZES:
# --- EXP-2.2, OPAQUE SNOW AND OPAQUE ICE
ax=axs[1,0]
h_biophytolte4_5_EXP2_2=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP2_2,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP2_2=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP2_2,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP2_2=ax.plot(
first_year,
array1d_iT_biophyto10_EXP2_2,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP2_2=ax.plot(
first_year,
array1d_iT_biophyto15_EXP2_2,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP2_2=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP2_2,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' d ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (E) EXP-2: SUM C_PHYTO: SIZES:
# --- EXP-2.1, OPAQUE SNOW AND TRANSLUCENT ICE
ax=axs[1,1]
h_biophytolte4_5_EXP2_1=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP2_1,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP2_1=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP2_1,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP2_1=ax.plot(
first_year,
array1d_iT_biophyto10_EXP2_1,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP2_1=ax.plot(
first_year,
array1d_iT_biophyto15_EXP2_1,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP2_1=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP2_1,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' e ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (F) EXP-2: SUM C_PHYTO: SIZES:
# --- EXP-0, TRANSLUCENT SNOW AND TRANSLUCENT ICE
ax=axs[1,2]
h_biophytolte4_5_EXP0=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP0,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP0=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP0,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP0=ax.plot(
first_year,
array1d_iT_biophyto10_EXP0,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP0=ax.plot(
first_year,
array1d_iT_biophyto15_EXP0,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP0=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP0,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' f ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- ADDITIONAL X-AXIS
ax=axs[2,0]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,1]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,2]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0,0].set_position([0.074,0.460, 0.292,0.372])
axs[0,1].set_position([0.388,0.460, 0.292,0.372])
axs[0,2].set_position([0.702,0.460, 0.292,0.372])
axs[1,0].set_position([0.074,0.066, 0.292,0.372])
axs[1,1].set_position([0.388,0.066, 0.292,0.372])
axs[1,2].set_position([0.702,0.066, 0.292,0.372])
axs[2,0].set_position([0.074,0.037, 0.292,0.020])
axs[2,1].set_position([0.388,0.037, 0.292,0.020])
axs[2,2].set_position([0.702,0.037, 0.292,0.020])
with plt.style.context('mplstyles/exp1_groups.mplstyle'):
# --- PLOT
fig,axs=plt.subplots(nrows=3,ncols=3,figsize=(29.445,12.825))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_exp2_groups.'+ext
plt.savefig(outfile,dpi=dpi)
Figure <sup_exp2_groups>. Sensitivity simulations by phytoplankton groups and size classes: EXP-2 with differing light under sea ice. Left column: EXP-2.2 with no light under snow-covered sea ice and no light under bare sea ice. Middle column: EXP-2.1 with no light under snow-covered sea ice but light under bare sea ice. Right column: EXP-0 with light under snow-covered sea ice and under bare sea ice. First row: Vertically integrated biomass for each functional group (0-100m). pico: picophytoplankton. other nano: other nanophytoplankton. dino: mixotrophic dinoflagellate. Second row: Vertically integrated biomass (0-100m) for different size classes or intervals of size classes. The size classes were described with their equivalent spherical diameter (ESD). Note that the interval of size classes <= 4.5 μm gathered together four size classes of the model and that the interval of size classes >= 20 μm gathered together eight size classes of the model. The vertical dotted line is the date of the snow melt completion. The grey rectangles represent the time of year before the sea ice break-up and after the sea ice freeze-up.
EXP-3¶
Differing pre-bloom nitrate.
EXP-3: infiles¶
infile_car_EXP3_times0_50\
='data/DataS6_output_mitgcm/exp3_times0_50/car.0000000000.t001.nc'
infile_car_EXP3_times1_00\
='data/DataS6_output_mitgcm/exp3_times1_00/car.0000000000.t001.nc'
infile_car_EXP3_times2_00\
='data/DataS6_output_mitgcm/exp3_times2_00/car.0000000000.t001.nc'
EXP-3: $\mathrm{ \sum C_{phyto} }$¶
$\mathrm{ \sum C_{phyto} }$ is the vertically integrated $\mathrm{ C_{phyto} }$ (0-100m) in $\mathrm{ mg\ C\ m^{-2} }$.
depth_end=100 # m
EXP-3: $\mathrm{ \sum C_{phyto} }$: groups¶
$\mathrm{ \sum C_{phyto} }$ for each phytoplankton group.
EXP-3: $\mathrm{ \sum C_{phyto} }$: groups: EXP-3, $\mathrm{ [NO_3]^{winter}*1/2 }$¶
array1d_iT_biodiatoms_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP3_times0_50\
=array1d_iT_biodiatoms_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biodino_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP3_times0_50\
=array1d_iT_biodino_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_bioother_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP3_times0_50\
=array1d_iT_bioother_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biopico_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP3_times0_50\
=array1d_iT_biopico_EXP3_times0_50/1000 # mg C -> g C
EXP-3: $\mathrm{ \sum C_{phyto} }$: groups: EXP-3, reference¶
array1d_iT_biodiatoms_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP3_times1_00\
=array1d_iT_biodiatoms_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biodino_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP3_times1_00\
=array1d_iT_biodino_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_bioother_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP3_times1_00\
=array1d_iT_bioother_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biopico_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP3_times1_00\
=array1d_iT_biopico_EXP3_times1_00/1000 # mg C -> g C
EXP-3: $\mathrm{ \sum C_{phyto} }$: groups: EXP-3, $\mathrm{ [NO_3]^{winter}*2 }$¶
array1d_iT_biodiatoms_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biodiatoms(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodiatoms_EXP3_times2_00\
=array1d_iT_biodiatoms_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biodino_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biodino(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biodino_EXP3_times2_00\
=array1d_iT_biodino_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_bioother_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_bioother(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_bioother_EXP3_times2_00\
=array1d_iT_bioother_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biopico_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biopico(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biopico_EXP3_times2_00\
=array1d_iT_biopico_EXP3_times2_00/1000 # mg C -> g C
EXP-3: $\mathrm{ \sum C_{phyto} }$: sizes¶
$\mathrm{ \sum C_{phyto} }$ for each size class.
EXP-3: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-3, $\mathrm{ [NO_3]^{winter}*1/2 }$¶
array1d_iT_biophytolte4_5_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP3_times0_50\
=array1d_iT_biophytolte4_5_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP3_times0_50\
=array1d_iT_biophyto6_6_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biophyto10_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP3_times0_50\
=array1d_iT_biophyto10_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biophyto15_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP3_times0_50\
=array1d_iT_biophyto15_EXP3_times0_50/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP3_times0_50\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP3_times0_50,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP3_times0_50\
=array1d_iT_biophytoge20_EXP3_times0_50/1000 # mg C -> g C
EXP-3: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-3, reference¶
array1d_iT_biophytolte4_5_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP3_times1_00\
=array1d_iT_biophytolte4_5_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP3_times1_00\
=array1d_iT_biophyto6_6_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biophyto10_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP3_times1_00\
=array1d_iT_biophyto10_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biophyto15_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP3_times1_00\
=array1d_iT_biophyto15_EXP3_times1_00/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP3_times1_00\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP3_times1_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP3_times1_00\
=array1d_iT_biophytoge20_EXP3_times1_00/1000 # mg C -> g C
EXP-3: $\mathrm{ \sum C_{phyto} }$: sizes: EXP-3, $\mathrm{ [NO_3]^{winter}*2 }$¶
array1d_iT_biophytolte4_5_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biophytolte4_5(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytolte4_5_EXP3_times2_00\
=array1d_iT_biophytolte4_5_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biophyto6_6_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biophyto6_6(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto6_6_EXP3_times2_00\
=array1d_iT_biophyto6_6_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biophyto10_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biophyto10(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto10_EXP3_times2_00\
=array1d_iT_biophyto10_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biophyto15_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biophyto15(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophyto15_EXP3_times2_00\
=array1d_iT_biophyto15_EXP3_times2_00/1000 # mg C -> g C
array1d_iT_biophytoge20_EXP3_times2_00\
=read_mitgcm.get_array1d_iT_biophytoge20(
infile=infile_car_EXP3_times2_00,
array1d_iT1y_iT=array1d_iT1y_iT,
array1d_idepth_delR=drF,
depth_end=depth_end
)
array1d_iT_biophytoge20_EXP3_times2_00\
=array1d_iT_biophytoge20_EXP3_times2_00/1000 # mg C -> g C
def make_plots(axs):
box=dict(pad=5,alpha=0)
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
locs=np.array([151,
181, 212, 243, 273, 304, 334])
labels=('Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
xmin=135 # 15 May (0-based)
xmax=364 # 31 December (0-based)
xlims=(xmin,xmax)
ymin=0
ymax=5
ylims=(ymin,ymax)
# --- (A) EXP-3: SUM C_PHYTO: GROUPS:
# --- EXP-3, [NO_3]^WINTER*1/2
ax=axs[0,0]
h_biodiatoms_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP3_times0_50,
'-',color='red',label='diatom')
h_biodino_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biodino_EXP3_times0_50,
'-',color='magenta',label='dino')
h_bioother_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_bioother_EXP3_times0_50,
'-',color='blue',label='other')
h_biopico_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biopico_EXP3_times0_50,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-3, $\mathrm{ [NO_3]^{winter}*1/2 }$',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (B) EXP-3: SUM C_PHYTO: GROUPS:
# --- EXP-3, REFERENCE
ax=axs[0,1]
h_biodiatoms_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP3_times1_00,
'-',color='red',label='diatom')
h_biodino_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biodino_EXP3_times1_00,
'-',color='magenta',label='dino')
h_bioother_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_bioother_EXP3_times1_00,
'-',color='blue',label='other')
h_biopico_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biopico_EXP3_times1_00,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-3, REFERENCE',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (C) EXP-3: SUM C_PHYTO: GROUPS:
# --- EXP-3, [NO_3]^WINTER*2
ax=axs[0,2]
h_biodiatoms_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biodiatoms_EXP3_times2_00,
'-',color='red',label='diatom')
h_biodino_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biodino_EXP3_times2_00,
'-',color='magenta',label='dino')
h_bioother_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_bioother_EXP3_times2_00,
'-',color='blue',label='other')
h_biopico_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biopico_EXP3_times2_00,
'--',color='green',label='pico')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' c ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
plt.text(0.5,1.2,'EXP-3, $\mathrm{ [NO_3]^{winter}*2 }$',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (D) EXP-3: SUM C_PHYTO: SIZES:
# --- EXP-3, [NO_3]^WINTER*1/2
ax=axs[1,0]
h_biophytolte4_5_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP3_times0_50,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP3_times0_50,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biophyto10_EXP3_times0_50,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biophyto15_EXP3_times0_50,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP3_times0_50=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP3_times0_50,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylabel('$\mathrm{ \sum\ C_{phyto} }$\n'
'($\mathrm{ g\ C\ m^{-2} }$)')
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' d ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (E) EXP-3: SUM C_PHYTO: SIZES:
# --- EXP-3, REFERENCE
ax=axs[1,1]
h_biophytolte4_5_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP3_times1_00,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP3_times1_00,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biophyto10_EXP3_times1_00,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biophyto15_EXP3_times1_00,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP3_times1_00=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP3_times1_00,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' e ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- (F) EXP-3: SUM C_PHYTO: SIZES:
# --- EXP-3, [NO_3]^WINTER*2
ax=axs[1,2]
h_biophytolte4_5_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biophytolte4_5_EXP3_times2_00,
':',color='cyan',label='$\mathrm{ \leq 4.5\ \mu m }$')
h_biophyto6_6_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biophyto6_6_EXP3_times2_00,
'--',color='orange',label='$\mathrm{ 6.6\ \mu m }$')
h_biophyto10_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biophyto10_EXP3_times2_00,
'-',color='red',label='$\mathrm{ 10\ \mu m }$')
h_biophyto15_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biophyto15_EXP3_times2_00,
'-',color='magenta',label='$\mathrm{ 15\ \mu m }$')
h_biophytoge20_EXP3_times2_00=ax.plot(
first_year,
array1d_iT_biophytoge20_EXP3_times2_00,
'-',color='black',label='$\mathrm{ 20\ \mu m }$')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.yaxis.set_ticks_position('both')
ax.set_ylim(ylims)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.05,0.90,' f ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
ax.axvline(obs_nosnow,color='black',linestyle=':')
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup-xmin,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=xmax-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25)
ax.add_patch(rect)
# --- ADDITIONAL X-AXIS
ax=axs[2,0]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,1]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- ADDITIONAL X-AXIS
ax=axs[2,2]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0,0].set_position([0.074,0.460, 0.292,0.372])
axs[0,1].set_position([0.388,0.460, 0.292,0.372])
axs[0,2].set_position([0.702,0.460, 0.292,0.372])
axs[1,0].set_position([0.074,0.066, 0.292,0.372])
axs[1,1].set_position([0.388,0.066, 0.292,0.372])
axs[1,2].set_position([0.702,0.066, 0.292,0.372])
axs[2,0].set_position([0.074,0.037, 0.292,0.020])
axs[2,1].set_position([0.388,0.037, 0.292,0.020])
axs[2,2].set_position([0.702,0.037, 0.292,0.020])
with plt.style.context('mplstyles/exp1_groups.mplstyle'):
# --- PLOT
fig,axs=plt.subplots(nrows=3,ncols=3,figsize=(29.445,12.825))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_exp3_groups.'+ext
plt.savefig(outfile,dpi=dpi)
Figure <sup_exp3_groups>. Sensitivity simulations by phytoplankton groups and size classes: EXP-3 prescribed pre-bloom nitrate concentrations. The nitrate concentrations before the bloom initiation ($\mathrm{ [NO_3]^{winter} }$) were prescribed at 1/2 and 2 times that in the reference simulation. First row: Vertically integrated biomass for each functional group (0-100m). pico: picophytoplankton. other nano: other nanophytoplankton. dino: mixotrophic dinoflagellate. Second row: Vertically integrated biomass (0-100m) for different size classes or intervals of size classes. The size classes were described with their equivalent spherical diameter (ESD). Note that the interval of size classes <= 4.5 μm gathered together four size classes of the model and that the interval of size classes >= 20 μm gathered together eight size classes of the model. The vertical dotted line is the date of the snow melt completion. The grey rectangles represent the time of year before the sea ice break-up and after the sea ice freeze-up.
Plot: mean Chl $a$ and $\mathrm{ C_{phyto} }$¶
def make_plots(axs):
box=dict(pad=5,alpha=0)
locs=np.array([0, 31, 59, 90, 120, 151,
181, 212, 243, 273, 304, 334])
labels=('Jan-1','Feb-1','Mar-1','Apr-1','May-1','Jun-1',
'Jul-1','Aug-1','Sep-1','Oct-1','Nov-1','Dec-1')
xlims=(first_year[0],first_year[-1])
xmin=xlims[0]
xmax=xlims[1]
# --- SNOW AND ICE THICKNESS
ax=axs[0]
ymin=-1.7
ymax=1.0
xsnow=(snow_thickness_df['doy']-1).to_numpy()
ysnow=snow_thickness_df['sample_thickness_m_average'].to_numpy()
h1=ax.plot(xsnow,ysnow,
'o',color='black',markersize=6,label='obs')
xsnow_to_fill=np.arange(xsnow[0],xsnow[-1]+1)
ysnow_to_fill=np.interp(xsnow_to_fill,xsnow,ysnow)
h2=ax.fill_between(xsnow_to_fill,ysnow_to_fill,
linestyle='-',color='aliceblue')
ax.text(95,0.3,'snow',
horizontalalignment = 'left',verticalalignment = 'center')
xice=(ice_thickness_df['doy']-1).to_numpy()
yice=-ice_thickness_df['sample_thickness_m_average'].to_numpy()
h3=ax.plot(xice,yice,'o',color='black',markersize=6)
xice_to_fill=np.arange(xice[0],xice[-1]+1)
yice_to_fill=np.interp(xice_to_fill,xice,yice)
h4=ax.fill_between(xice_to_fill,yice_to_fill,
linestyle='-',color='dodgerblue')
ax.text(95,-0.40,'ice',
horizontalalignment = 'left',verticalalignment = 'center')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylim(ymin,ymax)
ax.yaxis.set_ticks_position('both')
ax.set_ylabel('Thickness\n(m)',bbox=box)
ax.grid()
ax.legend(loc='upper right')
plt.text(0.03,0.75,' a ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=25)
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=364-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
ax.axvline(obs_nosnow,color='black',linestyle=':')
# --- MEAN CHL A AND C_PHYTO 0-40M
ax=axs[1]
ymin=7E-3
ymax=5E0
# mg Chl m^-3
h1=ax.plot(first_year,array1d_iT_modmeanchl,
'-',lw=3,color='black',label='model Chl')
h2=ax.plot(first_year,array1d_iT_obsmeanchlHPLC,
'o',markersize=6,color='black',label='obs Chl')
ax.set_yscale('log')
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels([])
ax.set_ylim(ymin,ymax)
ax.set_ylabel('Mean Chl $a$\n'
'(0-40m, $\mathrm{ mg\ Chl\ m^{-3} }$)')
ax.grid()
axt = ax.twinx()
# mg C m^-3
h3=axt.plot(first_year,array1d_iT_modmeancphyto,
'--',lw=3,color='red',label='model C')
axt.set_yscale('log')
axt.set_ylim(*[100*l for l in ax.get_ylim()])
axt.set_ylabel('Mean phyto biomass\n'
'(0-40m, $\mathrm{ mg\ C\ m^{-3} }$)')
# ref: https://stackoverflow.com/questions/5484922/secondary-axis-with-twinx-how-to-add-to-legend
# answer by Paul
# lns=h1+h2+h3
lns=h3+h1+h2
labs = [l.get_label() for l in lns]
ax.legend(lns,labs)
plt.text(0.03,0.93,' b ',
horizontalalignment = 'center',verticalalignment='center',
transform=ax.transAxes, fontsize=30)
rect=mpl.patches.Rectangle(xy=(xmin,ymin),
width=obs_breakup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
rect=mpl.patches.Rectangle(xy=(sim_freezeup,ymin),
width=364-sim_freezeup,height=ymax-ymin,
linewidth=0,facecolor='grey',alpha=0.25,
zorder=0)
ax.add_patch(rect)
ax.axvline(obs_nosnow,color='black',linestyle=':')
# --- ADDITIONAL X-AXIS
ax=axs[2]
ax.set_xlim(xlims)
ax.set_xticks(locs)
ax.set_xticklabels(labels)
ax.axes.get_yaxis().set_visible(False)
ax.spines['top'].set_visible(False)
# --- POSITION
axs[0].set_position( [0.12,0.75 ,0.75,0.179])
axs[1].set_position( [0.12,0.107 ,0.75,0.607])
axs[2].set_position( [0.12,0.061 ,0.75,0.032])
with plt.style.context('mplstyles/Chla_C_mean40m.mplstyle'):
# Plot
fig,axs=plt.subplots(3,1,figsize=(18.1,7.858))
make_plots(axs)
fig.align_ylabels(axs[:])
# --- SAVE
outfile='figures_timing_sup_progress/sup_Chla_C_mean40m.'+ext
plt.savefig(outfile,dpi=dpi)
Figure S<sup_Chla_C_mean40m>. Mean phytoplankton biomass and chlorophyll a (0-40m). (a) Observed snow and ice thickness. (b) Mean (0-40m) Chl a and biomass. Dots are in situ Chl a. The black solid line is the model Chl a. The red dashed line is the model phytoplankton biomass. The black vertical dotted line is the date of snow melt completion. The grey shading represents the time of year before the sea ice break-up and after the sea ice freeze-up.
Maps¶
Maps: Function¶
'''latlon2xy convenience function to convert latlon to chosen projection,
and figure object'''
def latlon2xy(lon,lat,proj):
xyz = proj.transform_points(ccrs.PlateCarree(),lon,lat).squeeze()
x,y,_ = np.split(xyz,3,axis=len(xyz.shape)-1)
x[np.isinf(x)]=np.nan
y[np.isinf(y)]=np.nan
return x.squeeze(),y.squeeze()
Maps: stations¶
# STATIONS
Station=collections.namedtuple(
'Station',
'lon lat marker markeredgecolor markerfacecolor markersize')
# Stations for the map of CCGS Amundsen 2018, leg 2b
station1=Station(-60.3919999999999,68.3096,'o','none','red',10)
station2=Station(-61.463,68.0358,'o','none','red',10)
station3=Station(-62.355,67.8668,'o','none','red',10)
station4=Station(-63.5639999999999,67.54,'o','none','red',10)
# Qikiqtarjuaq SIC
station5=Station(-63.783,67.482,'*','black','yellow',30)
station6=Station(-64.638,67.239,'o','none','red',10)
stations_am=[station1,station2,station3,station4,station5,station6]
# Stations for the map of Baffin Bay
station1=Station(station1.lon,station1.lat,station1.marker,
station1.markeredgecolor,station1.markerfacecolor,5)
station2=Station(station2.lon,station2.lat,station2.marker,
station2.markeredgecolor,station2.markerfacecolor,5)
station3=Station(station3.lon,station3.lat,station3.marker,
station3.markeredgecolor,station3.markerfacecolor,5)
station4=Station(station4.lon,station4.lat,station4.marker,
station4.markeredgecolor,station4.markerfacecolor,5)
# Qikiqtarjuaq SIC
station5=Station(station5.lon,station5.lat,station5.marker,
station5.markeredgecolor,station5.markerfacecolor,10)
station6=Station(station6.lon,station6.lat,station6.marker,
station6.markeredgecolor,station6.markerfacecolor,5)
stations_bb=[station1,station2,station3,station4,station5,station6]
Maps: Plot¶
plotbathymetry=True
with plt.style.context('mplstyles/maps.mplstyle'):
fig=plt.figure(figsize=(10, 16))
# --- CCGS AMUNDSEN 2018, LEG 2B
proj_am=ccrs.LambertConformal(central_longitude=-63,
central_latitude=67.5,
standard_parallels=(67,68))
ax_am=fig.add_subplot(223,projection = proj_am)
ax_am.coastlines(resolution='10m')
ax_am.add_feature(cartopy.feature.LAND)
if not plotbathymetry:
ax_am.add_feature(cartopy.feature.OCEAN)
MapExtent = [-66,-60,66,69]
ax_am.set_extent(MapExtent, ccrs.PlateCarree())
gl = ax_am.gridlines(crs=ccrs.PlateCarree(), draw_labels=False,
linewidth=1, color='gray', alpha=0.5,
linestyle='--')
gl.xlocator = mticker.FixedLocator(np.arange(-68,-58))
gl.ylocator = mticker.FixedLocator(np.arange(65,70))
# stations
for station in stations_am:
xstation,ystation = latlon2xy(np.array([station.lon]),
np.array([station.lat]),
proj_am)
marker=station.marker
markerfacecolor=station.markerfacecolor
markersize=station.markersize
markeredgecolor=station.markeredgecolor
ax_am.plot(xstation,ystation,marker,markeredgecolor=markeredgecolor,
markerfacecolor=markerfacecolor,markersize=markersize)
# longitude labels
y0,_ = ax_am.get_ylim()
for lon in [-66,-64,-62,-60]:
x0,_ = proj_am.transform_point(lon,MapExtent[2]-1,
src_crs = ccrs.PlateCarree())
plt.text(x0,y0-2e4,'{:2d}$^\circ$W'.format(-lon),
horizontalalignment = 'center',verticalalignment='top')
# latitude labels
x0,_ = ax_am.get_xlim()
someLons = np.arange(-80,-50)
for lat in np.arange(66,69+1):
# interpolate latitude circle to map boundary
xyzArray = proj_am.transform_points(ccrs.PlateCarree(),
someLons,lat*np.ones_like(someLons))
x = xyzArray[:,0]
y = xyzArray[:,1]
y0 = np.interp(x0,x,y)
plt.text(x0-2e4,y0,'{:2d}$^\circ$N'.format(lat),
horizontalalignment = 'right',verticalalignment='center')
if plotbathymetry:
# color bar for bathymetry
cmap = plt.cm.get_cmap('Greys')
cmap = cmap.from_list('Custom cmap',
[cmap(i) for i in range(150,50,-1)], 100)
cmap.set_over('white')
bounds = [0,500,1000,2000]
norm = mpl.colors.BoundaryNorm(bounds, cmap.N)
# IBCAO bathymetry
if plotbathymetry:
top = scipy.io.loadmat('data/DataS1_observations_IBCAO_1min_bathy.mat')
lon_ibcao,lat_ibcao=np.meshgrid(top['x'],top['y'])
z_ibcao=-top['IBCAO_1min']
hb_am=ax_am.contourf(*latlon2xy(lon_ibcao,lat_ibcao,proj_am),
z_ibcao,levels=[0,500,1000,2000],
cmap=cmap,norm=norm,extend='max')
axcb_am= fig.colorbar(hb_am,ax=ax_am,orientation='horizontal')
axcb_am.ax.set_title('Water depth (m)')
axcb_am.ax.set_position([0.70,-0.10, 0.25,0.20])
# labels for stations
for i,station in enumerate(stations_am):
xstation,ystation = latlon2xy(np.array([station.lon]),
np.array([station.lat]),
proj_am)
ax_am.text(xstation+3e3,ystation+1e4,i+1,
horizontalalignment = 'center',
verticalalignment='center')
# text
b_lon,b_lat = -66.5, 69.1,
ax_am.text(b_lon,b_lat,' b ',
horizontalalignment = 'left',fontsize=30,
transform=ccrs.PlateCarree())
# --- BAFFIN BAY
proj_bb=ccrs.LambertConformal(central_longitude=-60,
central_latitude=70,
standard_parallels=(65,75))
ax_bb=fig.add_subplot(221,projection=proj_bb)
ax_bb.coastlines(resolution='50m')
ax_bb.add_feature(cartopy.feature.LAND)
if not plotbathymetry:
ax_bb.add_feature(cartopy.feature.OCEAN)
MapExtent = [-76,-48,62,79]
ax_bb.set_extent(MapExtent, ccrs.PlateCarree())
gl = ax_bb.gridlines(crs=ccrs.PlateCarree(), draw_labels=False,
linewidth=1, color='gray', alpha=0.5,
linestyle='--')
gl.xlocator = mticker.FixedLocator([-80,-75,-70,-65,-60,
-55,-50,-45,-40,-35])
gl.ylocator = mticker.FixedLocator(np.arange(58,84,3))
# stations
for station in stations_bb:
xstation,ystation = latlon2xy(np.array([station.lon]),
np.array([station.lat]),
proj_bb)
marker=station.marker
markerfacecolor=station.markerfacecolor
markersize=station.markersize
markeredgecolor=station.markeredgecolor
ax_bb.plot(xstation,ystation,marker,markeredgecolor=markeredgecolor,
markerfacecolor=markerfacecolor,markersize=markersize)
# longitude labels
y0,_ = ax_bb.get_ylim()
for lon in [-70,-60,-50]:
x0,_ = proj_bb.transform_point(lon,MapExtent[2]-1,
src_crs = ccrs.PlateCarree())
plt.text(x0,y0-2e4,'{:2d}$^\circ$W'.format(-lon),
horizontalalignment = 'center',verticalalignment='top',
fontsize=20)
# latitude labels
x0,_ = ax_bb.get_xlim()
someLons = np.arange(-120,-40,1)
for lat in np.arange(64,78,3):
# interpolate latitude circle to map boundary
xyzArray = proj_bb.transform_points(ccrs.PlateCarree(),
someLons,lat*np.ones_like(someLons))
x = xyzArray[:,0]
y = xyzArray[:,1]
y0 = np.interp(x0,x,y)
plt.text(x0-2e4,y0,'{:2d}$^\circ$N'.format(lat),
horizontalalignment = 'right',verticalalignment='center')
# IBCAO bathymetry
if plotbathymetry:
hb_bb=ax_bb.contourf(*latlon2xy(lon_ibcao,lat_ibcao,proj_bb),
z_ibcao,levels=[0,500,1000,2000],
cmap=cmap,norm=norm,extend='max')
# text
a_lon,a_lat = -98, 77.5
ax_bb.text(a_lon,a_lat,' a ',
horizontalalignment = 'left',fontsize=30,
transform=ccrs.PlateCarree())
greenland_lon,greenland_lat = -38, 77
ax_bb.text(greenland_lon,greenland_lat,'Greenland',
horizontalalignment = 'right',
transform=ccrs.PlateCarree())
baffinbay_lon,baffinbay_lat = -63, 72
ax_bb.text(baffinbay_lon,baffinbay_lat,'Baffin\nBay',
horizontalalignment = 'right',
transform=ccrs.PlateCarree())
# rectangle
extentStudyarea = ax_am.get_extent()
vertices = np.array([
[extentStudyarea[0], extentStudyarea[2]], # left, bottom
[extentStudyarea[0], extentStudyarea[3]], # left, top
[extentStudyarea[1], extentStudyarea[3]], # right, top
[extentStudyarea[1], extentStudyarea[2]], # right, bottom
[extentStudyarea[0], extentStudyarea[2]]
])
xVert,yVert = vertices.transpose()
xBox = np.interp(np.arange(0,4,.02),np.arange(0,4.1,1),xVert)
yBox = np.interp(np.arange(0,4,.02),np.arange(0,4.1,1),yVert)
xyzBoxLatlon = ccrs.PlateCarree().transform_points(ax_am.projection,
xBox,yBox)
xBoxLatlon = xyzBoxLatlon[:,0]
yBoxLatlon = xyzBoxLatlon[:,1]
xyzBoxTransformed = proj_bb.transform_points(ccrs.PlateCarree(),
xBoxLatlon,yBoxLatlon)
xBoxTransformed = xyzBoxTransformed[:,0]
yBoxTransformed = xyzBoxTransformed[:,1]
ax_bb.plot(xBoxTransformed,yBoxTransformed,'r')
ax_bb.relim()
# --- WORLD MAP
proj_wm = ccrs.NearsidePerspective(central_longitude=-60,
central_latitude=60)
ax_wm = fig.add_subplot(222,projection = proj_wm)
ax_wm.set_global()
ax_wm.add_feature(cartopy.feature.LAND, zorder=0, edgecolor='black')
ax_wm.add_feature(cartopy.feature.OCEAN)
gl = ax_wm.gridlines(crs=ccrs.PlateCarree(), draw_labels=False,
linewidth=1, color='gray', alpha=0.5,
linestyle='--')
extentStudyarea = ax_bb.get_extent()
vertices = np.array([
[extentStudyarea[0], extentStudyarea[2]], # left, bottom
[extentStudyarea[0], extentStudyarea[3]], # left, top
[extentStudyarea[1], extentStudyarea[3]], # right, top
[extentStudyarea[1], extentStudyarea[2]], # right, bottom
[extentStudyarea[0], extentStudyarea[2]]
])
xVert,yVert = vertices.transpose()
xBox = np.interp(np.arange(0,4,.02),np.arange(0,4.1,1),xVert)
yBox = np.interp(np.arange(0,4,.02),np.arange(0,4.1,1),yVert)
xyzBoxLatlon = ccrs.PlateCarree().transform_points(ax_bb.projection,
xBox,yBox)
xBoxLatlon = xyzBoxLatlon[:,0]
yBoxLatlon = xyzBoxLatlon[:,1]
xyzBoxTransformed = proj_wm.transform_points(ccrs.PlateCarree(),
xBoxLatlon,yBoxLatlon)
xBoxTransformed = xyzBoxTransformed[:,0]
yBoxTransformed = xyzBoxTransformed[:,1]
ax_wm.plot(xBoxTransformed,yBoxTransformed,'r')
ax_wm.set_position( [0.60,0.60 ,0.38,0.38])
ax_bb.set_position( [0.10,0.53 ,0.60,0.42])
ax_am.set_position( [0.10,0.05 ,0.60,0.42])
# --- SAVE
if plotbathymetry:
plt.savefig('figures_timing_sup_progress/sup_maps.'+ext,dpi=dpi)
else:
plt.savefig('figures_timing_sup_progress/sup_maps_nobathy.'+ext,dpi=dpi)
/var/folders/vq/5y09539n705gpb79d6ysx09r0000gr/T/ipykernel_78452/1440760513.py:57: MatplotlibDeprecationWarning: The get_cmap function was deprecated in Matplotlib 3.7 and will be removed two minor releases later. Use ``matplotlib.colormaps[name]`` or ``matplotlib.colormaps.get_cmap(obj)`` instead. /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/constructive.py:181: RuntimeWarning: invalid value encountered in buffer /Users/maximebenoit-gagne/anaconda3/envs/timing/lib/python3.8/site-packages/shapely/predicates.py:798: RuntimeWarning: invalid value encountered in intersects
Fig. <sup_maps>. Maps of Baffin Bay and the study site. The yellow star marker represents the Qikiqtarjuaq sea ice camp location (Oziel et al., 2019). (a) Map of Baffin Bay. (b) Map of the six stations of the leg 1b of the CCGS Amundsen expedition in 2018 including the location of the Qikiqtarjuaq ice camp.
Plot: T-S diagram¶
Data is from CCGS Amundsen 2018, leg 2b. The location is shown on Fig. <sup_maps> (maps_bathymetry.jpeg). Code by Raphaël Larouche and Maxime Benoît-Gagné.
def density_contour(ax, sgrid, tgrid):
"""
Plot of density contour.
:param ax:
:param sgrid:
:param tgrid:
:return:
"""
dens = gsw.sigma0(*np.meshgrid(sgrid, tgrid))
c_dens = ax.contour(sgrid, tgrid, dens, levels=np.arange(22, 29, 1), linestyles='dashed', colors='gray', linewidths=.8)
ax.clabel(c_dens, inline=True, fmt="%1.1f $\mathrm{kg\ m^{-3}}$", fontsize=10)
# Open data
path='data/DataS9_observations_CCGS_Amundsen/1802_ODV_truncated.csv'
data=pd.read_csv(path,header=[0],sep=',')
data=data.set_index("Station") # Set index of dataframe to station number
# Dictionary for station number
station_dic = {"station 1": 7,
"station 2": 9,
"station 3": 10,
"station 4": 12,
"station 5": 11,
"station 6": 13}
def make_plots(ax):
# Colormap for depth
colmap = mpl.cm.get_cmap("hot")
normalization = mpl.colors.Normalize(vmin=0.0, vmax=1100.0)
# Color each station
cl = {"station 1": "g",
"station 2": "b",
"station 3": "y",
"station 4": "c",
"station 5": "k",
"station 6": "m"}
# Linestyle each station
lstyle = {"station 1": "-",
"station 2": "-.",
"station 3": ":",
"station 4": "--",
"station 5": "-",
"station 6": "-."}
# Density lines
sgrid = np.linspace(29, 35, 100)
tgrid = np.linspace(-2, 2, 100)
density_contour(ax, sgrid, tgrid)
# Freeze curve
Sf = np.linspace(31, 34, 10)
Tf = gsw.CT_freezing(Sf, 0, 0)
ax.plot(Sf, Tf, "-", color="gray", linewidth=5)
# Loop
for key in station_dic.keys():
station_data=data.loc[station_dic[key]]
# Longitude, latitude
lon=station_data['Longitude [degrees_east]'].iloc[0]
lat=station_data['Latitude [degrees_north]'].iloc[0]
# Depth and normalization
depth=np.array(station_data["Depth [metres]"]).astype(list)
color_depth = [colmap(normalization(value)) for value in list(depth)]
# Salinity, conservative temperature data
salinity = np.array(station_data["Asal [g/kg]"]).astype(float)
cons_temp = np.array(station_data["Cont [deg C]"]).astype(float)
# Label
lab = "{0} ({1:.3f}$\degree$ W)".format(key, -1*lon)
# Plot data
if key == "station 5":
ax.plot(salinity, cons_temp, linewidth=5, color=cl[key], linestyle=lstyle[key], label=lab)
else:
ax.plot(salinity, cons_temp, linewidth=3, color=cl[key], linestyle=lstyle[key], label=lab)
# Figure TS diagram
ax.set_xlim(29, 35)
ax.set_ylim(-2, 2)
ax.set_yticks(np.arange(-1.5,2.5,0.5))
ax.set_ylabel('$\mathrm{ T_C\ (\degree C) }$')
ax.set_xlabel('$\mathrm{ S_A\ (g\ kg^{-1}) }$')
ax.legend(loc="best")
with plt.style.context('mplstyles/ts.mplstyle'):
# Plot
fig=plt.figure(figsize=(16, 16))
ax=fig.add_subplot(111)
make_plots(ax)
# --- SAVE
plt.savefig('figures_timing_sup_progress/sup_ts.'+ext,dpi=dpi)
/var/folders/vq/5y09539n705gpb79d6ysx09r0000gr/T/ipykernel_78452/2180134486.py:3: MatplotlibDeprecationWarning: The get_cmap function was deprecated in Matplotlib 3.7 and will be removed two minor releases later. Use ``matplotlib.colormaps[name]`` or ``matplotlib.colormaps.get_cmap(obj)`` instead.